arXiv: 1506.04565v2 [math.AP] 1 Dec 2015 


OPTIMAL TRANSPORTATION THEORY WITH REPULSIVE COSTS 


SIMONE DI MARINO, AUGUSTO GEROLIN, AND LUCA NENNA 


Contents 


Introduction ^ 

Organization of this Survey [2] 

What is new? [3| 

1. Why Multi-marginal Transport Theory for repulsive costs? [3| 

1.1. Brief introduction to Quantum Mechanics of A-body systems [3| 

1.2. Probabilistic Interpretation and Marginals [6| 

1.3. Density Functional Theory (DFT) [6| 

1.4. “Semi-classical limit” and Optimal Transport problem [8| 

2. DFT meets Optimal Transportation Theory fTOl 

2.1. Couplings and Multi-marginal Optimal Transportation problem 1101 

2.2. Multimarginal Optimal Transportation Problem 1101 

2.3. Dual formulation 1121 

2.4. Geometry of the Optimal Transport sets 1141 

2.5. Symmetric Case 1151 

3. Multi-marginal OT with Coulomb Cost 1161 

3.1. General theory: duality, equivalent formulations and many particles limit 1161 

3.2. The Monge problem: deterministic examples and counterexamples 1181 

4. Multi-marginal OT with repulsive Harmonic Cost 1191 

5. Multi-marginal OT for the Determinant 1241 

6. Numerics 117 ] 

6.1. The regularized problem and the Iterative Proportional Fitting Procedure [27| 

6.2. Numerical experiments: Coulomb cost 1301 

6.3. Numerical experiments: repulsive Harmonic cost 1301 

6.4. Numerical experiments: Determinant cost 1321 

7. Conclusion 1341 

Acknowledgements 1341 

References 1341 


Abstract. RICAM survey repulsive cost in multi-marginal optimal transportation theory 
and predictions for the future! 


Introduction 

The main objective of this survey is to present recent developments in optimal transport theory 
for repulsive costs with finitely many marginals. We are interested in characterizing, for a class 
of repulsive cost functions, the minimizers of the multi-marginal optimal transportation problem 

{MK.) mini f c(a:i,..., a:jv)d7 : 7 S ,...,/xjv) >, 

[ Jk<‘" J 

where c : ^ R U { 00 } is a lower semi-continuous cost function and n(R‘^^,/ri,... ,/rjv) 

denotes the set of couplings, i.e. probability measures ^{x \,... ,Xjv) having {li as marginals, for 
all 2 = 1,..., Y. c is, typically, a sum of convex decreasing functions, as the Coulomb interaction 

cost [niilTiiiiiiiiiEgiS] 

N N ^ 

c(xi,...,xjv) = ^ I _ 
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which is the case of most interest in DFT (see section [y. Also related to Quantum Mechanics, 
in the so-called weak interaction regime |12l 13811731 1851 T we will consider the repulsive harmonic 
potential as cost function c, 

N N 

c„(a;i,.. . ,xjv) = - ^ ^ \xj-Xi\^. 

i=lj=i+l 

In the last section, we focus on the multi-marginal problem for the determinant cost 
c{xi ,..., xn) = det(xi,..., Xiv), or c{xi ,. .., xn) = \ det(xi ,... ,X]\j)\, 
studied by Carlier and Nazaret in Ezi. 

Multi-marginal optimal transport showed to be powerful in many applications as, for instance, 
in economics |15l 1211 \22[ 1351171| . mathematical finance [HI 171 l33l 1341139| . image processing |87| . 
tomography ^ and statistics m- 

Although very natural in such applied problems, optimal transport problems with 3 or more 
marginals was only considered in literature during the 90’s by Knott and Smith, Olkin and Rachev 
and Riischendorf and Uckelmann [531 Ibfil (7^ . and later by Gangbo and Swiech |41| . In |41| . the 
authors proved the existence and uniqueness of Monge-type solution for the problem (Ad/C) for 
the quadratic cost c{xi ^... \^3 ~ 

A systematic study of the structure of transport plans and general conditions in order to 
guarantee the existence of Monge-type solutions for the problem (Ad/C) was made by Pass EZl 
El [69] and Kim and Pass [UlSo]. Several theoretical developments were made last years, as the 
study of partial multi-marginal transport [52], symmetric Monge-Kantorovich problem [HEllHll 
14311441 1611 l62l I89| and duality [ 6 | 1301 l48l I89| . Applications in Geometry and Analysis as the 

Wasserstein Barycenter problem and decoupling a system of PDE’s can be found in 1441 fST) . 

For a complete presentation of the multimarginal optimal transport theory with general costs 
functions, including a discussion of many applications, we refer to the survey ED- 

In this paper, we are rather analyzing the Multimarginal Optimal Transport Theory starting 
from the point of view of researchers working on Density Functional Theory. 

Paola Gori-Giorgi, Michel Seidl et al were interested in computing accurately the ground 
state of a molecule with a large number of atoms. As we will discuss later, after doing suitable 
approximations, the solutions of the optimal transport (Ad/C) for Coulomb cost (more precisely, 
the optimal potential), or more generally repulsive costs, can be seen as a good ansatz in order 
to construct solutions for Hohenberg-Kohn functional ( |1.12[ |. 

From a mathematical perspective the interesting case is when all marginals of the coupling 7 
in (Ad/C) are absolutely continuous with respect to Lebesgue measure and are all the same. In 
particular, since the cost functions are “repulsive”^ if Monge-type solutions exists, they should 
follow the rule “the further, the better!”, which means that we want to move the mass as much 
as we can. In other words, in the present case, optimal transport plans tend to be as spread as 
possible. But when a cost can be considered repulsive? It is difficult to find a good notion for 
this concept since usually it is based on the intuition that comes solving the particular problem. 
However a possible formal definition for repulsive cost is to require a mild condition, namely 
that we don’t want any two particles in the same position; it should be cheaper to have them in 
different position. For example can write the inequality for the particles 1 and 2: 

c(xi, xi, X 3 ,... ,a:iv) + c{x2,X2,xs,... ,xn) > 

c(xi,X2,Xs, ... ,Xn) c{x2,Xi,Xs, . . .,Xn)- 

Then, also all other inequalities for i ^ j should be considered. This condition is satisfied by all 
the considered costs, even the ones that doesn’t seem so “repulsive” at first, such as those of the 
form f{xi -h ... + Xn) with / convex function. 

Organization of this Survey. In section^ we present the physical motivation in considering 
repulsive costs and give an introduction of Density Functional Theory from a mathematical point 
of view. The main goal of this section is to explain how optimal solutions of Monge-Kantorovich 
problem can be useful to compute the ground state energy of a quantum system of N electrons. 

Section [^introduces briefly the general theory of multimarginal optimal transport, including a 
discussion about the Kantorovich duality and the notion of c-cyclical monotonicity applied to this 
framework. In addition, we give some highlights in the interesting case of symmetric transport 
plans. 




OPTIMAL TRANSPORTATION THEORY WITH REPULSIVE COSTS 


In sections and respectively, we discuss the recent development for the specific cases of 
transportation problem with Coulomb, repulsive harmonic and determinant cost. 

Finally, in section]^ we introduce an iterative method, recently used in many applications of 
optimal transport [9l|l0|, and we present some numerical experiments for the problems studied in 
previous sections. 

What is new? We present in section]^ some new results we obtained for cost functions 

c{xi,.. . ,xjs[) = f{xi + • • • + ^ R convex^ 

which, in particular, include the case of repulsive harmonic costs (see corollary |4.8| |. In addition, 
we show: 

(i) Existence of Fractal-like cyclical Monge solutions: In Theorem |4.6| we show the existence 

of couplings 7 = concentrated on the graphs of functions 

where /i a uniform probability measure on the d-dimensional cube. 
Moreover, T is not differentiable in any point, which differs substantially to the classical 
regularity theory in the classical 2 marginals case [131 [g. 

(ii) Optimal Transport plans supported on sets of Hausdorff dimension > 1; In standard 
optimal transportation (the 2 marginals case) the target is to understand if solutions of 
(AiK) are well-known to be Monge-type solutions: optimal transport plans are supported 
by the graph of a function T. In Section]^ we give a concrete example (see example |4.14[ l 
where optimal transport plans are supported on 2-dimensional Hausdorff sets. 

(iii) Numerical results: we present some interesting results for the repulsive harmonic and for 
the determinant cost. 

1. Why Multi-marginal Transport Theory for repulsive costs? 

1.1. Brief introduction to Quantum Mechanics of A^-body systems. Energies and 
geometries of a molecule depend sensibly on the kind of atom on a chemical environment. Semi- 
empirical models as Lewis structure explain some aspects, but they are far from being satisfactory, 
mostly because they are not quantitative. 

Classical behavior of atoms and molecules is described accurately, at least from a theoretical 
point of view, by quantum mechanics. However, in order to predict the chemical behavior of 
a molecule with a large number of electrons we need to deal with computational aspects and 
approximations are needed. For example, if we use a direct approximation of 10 grid points in 
R of the time-dependent Schrodinger equation in order to simulate the chemical behavior of the 
water molecule (H 2 O), which has 10 electrons and its position is represented in the space R^®, 
then we need to compute 10^® grid points. 

Quantum Chemistry studies the ground state of individual atoms and molecules, the excited 
states and the transition states that occur during chemical reactions. Many systems use the 
so-called Born-Oppenheimer approximation and many computations involve iterative and other 
approximation methods. 

The main goal of Quantum Chemistry is to describe accurately results for small molecular 
systems and to increase the size of large molecules that can be processed, which is very limited 
by scaling considerations. 

It means that from a Quantum Model, for which the only input is the atomic number of a 
molecule, we want to predict the evolution of a molecular system. Classical models consider a 
molecule with N electrons and M nuclei. We denote by a point Xi G R^ the position coordinates 
of the 2 -th electron of mass me and by Si its 2 -th electron spin. The charges, the masses and 
the positions of the o-th nucleus are represented, respectively, by G R-^^, E and 

Rc e (lR3)Af 

The state of the system is described by a time-dependent wave function i/i S T^([0,T];(S)7,^L^(R®x 
X tp{xi, si,..., Xff, sjv; -Rii ■ ■ • I Rm^R)- The group 1? = {tU} represents the 

spin of a particle. 

We say that a wave function is antisymmetric or fermionic if it changes sign under a simul¬ 
taneous exchange of two electrons i and j with the space coordinates Xi and Xj and spins Si and 
Sj, that is, 

• >a:^(Ar).«a(JV)) = sign(o-)V’(xi, si,. .. , Xjy, S]v) tr S Yjv 
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where yjsi denotes the permutation set of N elements. 

In physics literature, the Fermionic wave function obeys the Fermi-Dirac statistics and, in 
particular, the Pauli exclusion principle. On the other hand, a wave function is called symmetric 
or bosonic when it obeys Bose-Einstein statistics which also implies that when one swaps two 
bosons, the wave function of the system is unchanged. 

If it is not mentioned explicitly, we will suppose that the wave function is spinless, i.e. i/> is a 
function depending only on time, electrons and nuclei position , X 2 ,..., xn, Ri ,..., Rm)- 

The evolution of a wave function ijj boils down to a (time-dependent) Schrodinger Equation 

idt'lp = H'lp , H = Tn + Te + Vne + Vee + Vnn ( 1 * 1 ) 


where the operators Tn,Te, Vne, Vee and Vnn are defined by 
M 


Tn = —h - ^Ra (Euclei Kinetic Energy) 

Te = —hy -Aa;v (Electron Kinetic Energy) 

2me 


i=l 

N 


i=l Vq = 1 I** -““I 


i=l 

N N 

Vee = E E 

i=l j = l,j^i 
M M 


Potential Energy of Inte¬ 
raction Nuclei-Electron 


\Xi - Xj\ 


(Interaction Electron-Electron) 
(Interaction Nuclei-Nuclei) 


Born-Oppenheimer Approximation: The Born-Oppenheimer approximation (named for 
its original inventors. Max Born and Robert Oppenheimer) is based on the fact that nuclei are 
several thousand times heavier than electrons. The proton, itself, is approximately 2000 times more 
massive than an electron. Roughly speaking, in Born-Oppenheimer approximation we suppose 
that nuclei behave as point classical particles. This is reasonable, since typically 20007ne < ma < 
lOOOOOme. 

In other words, we suppose ^ rrie, we fix mg = 1 and we consider an ansatz of type 
'ilj{x^R,t) = 'ifeix, R)x{R,t), which physically means that the dynamics of the electron -ipe are 
decoupled of the dynamics of the atomic nuclei y. Substituting that ansatz in mi , we can 
show that the electronic part ipe{x, R) solves the following eigenvalue problem, so-called Eletronic 
Schrodinger Equation 

He'^e — AeV’ei He — Tg + Vne T Vee T Vnn (l*^) 

where the eigenvalue Ag = Ae(R) depends on the position vector R = (Ri,..., Rm) of th® atomic 
nuclei. The function He = He{R) is called Eletronic Hamiltonian. From the other side, through a 
formal argument, the nuclei wave function a solution of a Schrodinger Equation restricted 

to the nuclei with potential energy Ag, 

idtxiH) = {Tn + Xe(R))x(R)- 

We refer to m for all computations and formal deduction of those formulas. Another conse¬ 
quence of that approximation is that the nuclear components of the wave function are spatially 
more localized than the electronic component of the wave function. In the classical limit, the 
nuclei are fully localized at single points representing classical point particles. 


In Quantum Chemistry it is interesting to study the so called Geometric Optimization Problem. 
The Geometric Optimization Problem: Compute the following minimizer 

M M 


{IVl IVi 

E{Ri,..., Rm) + E E 


Za Z, 


ol^I3 


l/3=l,/37^a 1-^“ 


: (Ri,..., Rm) S 


j3M 


(1.3) 


where 


• The second term is the nuclei-nuclei iteration Vnn, as defined in the 






OPTIMAL TRANSPORTATION THEORY WITH REPULSIVE COSTS 


• The function E{R) = E{Ri, ..., Rm) corresponds to the effective potential created by 
the electrons and is, itself, given by a minimization problem (see Electronic minimization 
problem bellow). 

The value of the potential Eq = Eo{Ri,. .., Rm) for given fixed positions of the nuclei is 
obtained by solving the electronic problem: 

The Electronic minimization problem: Compute the lowest eingenvalue ( ^‘ground state 
energy”) Eq of the following linear operator, called Electronic Hamiltonian 

Definition 1.1 (Eletronic Hamiltonian). The Eletronic Hamiltonian is a linear operator H^i : 

^2 ^ ^ N N 

i = l i = l = 

where denote the set of square-integrable antisymmetric functions -0 : —>■ R, 

: R^ ^ R is a L^(R^-^) function and / : R —>■ R is a continuous function. 


The first term of the Eletronic Hamiltonian is the Kinetic Energy and the second term, the 
function vji(xi), depends only on the single electron. Typically, vji(xi) = — XlaLi —j- is 

the interaction electron-nuclei energy of the single electron. The function / depends only on 
the distance of two electrons and measures the electron-electron potential interaction; typically 
/ is the Coulomb Interaction f{xj — Xi) = l/\xj — x.i\ or the repulsive harmonic interaction 
f{xj — Xi) = —\xj — Xi\^. From both mathematical and physical viewpoint, it may be also useful 
to study more general convex and concave functions. 

The ground state of the functional EH is given by the Rayleigh-Ritz variational principle 

Eo = E[Ri,- ■ ■ ,Rm) = min{(V',-f/eii/'> ; V’ S H}, 

where "H = {ip S : ip is antisymmetric and Hi/jH = 1}. Equivalently, 

Eo = min|Te[tl'] + V„e[ip] + Veelip] ■■ Ip e "nj (1.4) 

where Te['0] is the Kinetic energy^ 

- f r 

n[ip\ = — T / •■■ T / YWxii’{xi,s\ ... ,xi^,sj^)\^dxi .. .dxj^ 


Vne^^] is the electron-nuclei interaction energy, 

N 


V„e[ip]= Y / ■■■ Y / ■ >a^iV.Siv)l^da;i 


... dx A 


and, Vee['0] is the electron-electron interaction energy 

N N 


Ve 


ee[iy] — 


E L - E /3E E 

i6z2-'k 


1 


-\ip{xi, SI , Xff, sjv)|^dxi .. . dx]^ 


Here, the ground state quantum A^-body problem refers to the problem of finding equilibrium 
states for system of type A theorem proved first by Zhislin [90], guarantees the existence of 

the minimum in ( |1.4| t. Some variants can be found in the literature due, for instance, to Lieb ^ 
Simon ESI. Lions I60I and Friesecke EZj. 

In computational chemistry, the Electronic minimization problem methods are used to compute 
the equilibrium configuration of molecules and solids. However, the most computationally practi¬ 
cable methods in this context are not numerical methods (in the sense in which this terminology 
is used in mathematics), but power series of analytical solutions of reduced models. 
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1.2. Probabilistic Interpretation and Marginals. The square norm of the wave function 
- ■ ■ ,xn,S];^) can be interpreted as a A^-point probability density distribution for the 
electrons are in the points Xi with spins 


/ 


\-ip{xi, Si, 


,3;jv,sjv)|^cix = 1 . 


(1.5) 

(1.6) 


If i/i is a L^(R^^) function we can define the single particle density by 

p{xi) = N 7]v(a:i, SI, ■ ■ • , xjy, S]v)da;2, • ■ ■ dx^, 

JrHn-1) 

where 7 jv represents the A^-point position density 

7jv(xi, • • ■ ,xjv) = ^ |i/)(a;i,si,--- ,xjv,Siv)p. 

Analogously, we can define the A)-density 

7fe(a:i,...,Xfe) = [ 7jv(a;i,si,--- ,a:jv,sjv)dxj,+i.---(i3:]v (1.7) 

We remark that when the particle is spinless the single particle density can be simply be 

written as 

p(xi) = N \ll){xi, XN)\'^dX2, ■ ■ ■ dxjq. 

7j13(N-1) 


The relevance of p and 72 to compute the ground state energy ( |1.4[ | is due to the fact that the 
interation nuclei-nuclei energy Vne['4’] and electron-electron energy \ 4 e[‘iA] depends only on these. 


Lemma 1.2 ([27|)' G x -tp antisymmetric and \\ip \\2 = 1? then 

Vne[ip] = [ v{x)p{x)dx, Vee[V’] = / 1—“— n2ix,y)dxdy. 

7 m 3 J^6 \x - y\ 


It is natural to wonder about the space of those densities arising from such a density -ip is 
defined, i.e, the space A of densities p : R^ —>■ R verifying for an antisymmetric ip such that 

ip,^ip G L^((R^ X Z^)^,R) and \\ip \\2 = 1- An explicit characterization exists (due to Lieb |57|h 


A = 


l^p : R^ —>■ R : p > 0, ^/p G and 



(1.8) 


1.3. Density Functional Theory (DFT). The Density Functional Theory is the standard 
approximation to quantum mechanics in simulations of a system with more than a dozen electrons 
and showed to be successful in many instances but has rare drastic failures as, for instance, in 
predicting the behavior of Cr 2 molecules m- DFT theory approximates quantum mechanics via 
variational principles for the marginal density 

p{xi) := / \^{xi,...,XM)\'^dx 2 ...dxM, (1.9) 

yi^d(iv-i) 

where ip is s. wave function associated to the W-body quantum problem as in equation 
Roughly speaking, DFT models are semi-empirical models of the pair density 

72(a;i,X2):= / \iIj{xi, ... ,XN)fdx3 .. .dx^- (1-10) 

yj{3(N-2) 

in terms of its marginal p. We simply write ip ^ p and ip 'y 2 to denote the relation between 
Ip and p, and ip and 72 . This means, respectively, that ip has single particle density p and pair 
density 72 . 

Concerning the ground state problem ( |1.4[ |, we notice that, after the lemma mi we are able 
to write the Electronic minimization problem l |1.4[ ) as 

Eq = inf |f^k[p]+ [ v{x)p(x)dx\, (l-H) 

I Jr3 ) 

with 

F"^[p] = inf|r[V>] + Vee[ 72 ] :pe (1-12) 
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where Fjjk is the so called Hohenherg-Kohn or Levy-Lieb functional. For a mathematical sake of 
completeness we also consider a weak repulsive interaction given by 


Ve4p] 





, XN)dx\ .. .dxN, 


(1.13) 


in a compact set K C and +00 otherwise. 

The next theorem is the main result in Density Functional Theory. 


Theorem 1.3 (Hohenberg-Kohn, Levy-Lieb). Let p E P(R^) be a probability density such that 
p E and / : R ^ R 6e a continuous function. There exists afunctional : 7^(R^) R 

depending only on the single-particle density p such that for any potential Vnm the exact Quantum 
Mechanics ground state energy satisfies 

_Eo = i^in I + W J vji(x)p{x)dx^ : p E T’(R^)| (1-14) 

Moreover, F^^[p] is given itself as a minimum problem 

= min ] (V", (- + /(Xj -Xi))'!/’) : 

where, ip ^ p means that ip has single-particle density p and A is the set 0 / p : R^ ^ R such that 
P ^ Oj y/p £ Lf^(R^) and p{x)dx = N. 

The Hohenberg-Kohn theorem states that the functional F^^ is universal in the sense that 
it does not depend on the molecular system under consideration. From a Physics point of view, it 
garantees that in a molecular system of N electrons, the single electron density p determines the 
pair density of the system (see corollary |1.4[ |. From a mathematical perspective, the proof of this 
theorem is a functional analysis exercise that, for sake of completeness, we are going to present in 
a version we learnt from Gero Friesecke. 


Proof. Firstly we remark that the non-universal part of the energy only depends on p: 

N N 

{'4),y^v(xi)ip)=j y'i;(xi)|V’(xi,... ,X]v)p(ix = Af / v{x)p(x)dx 

Then, we rewrite problem ( |1.4| l as a double inf problem 

Eo = inf pij), He; Ip) + N J v{r)p(r)dr'j 

= inf ^ inf pip,}lf,i ip)'j + N J v(r)p(r)di 
= inf + N J v{r)p{r)di 

where ((V’,He; V’>)- 


The present version stated in theorem |1.3| is due by Levy and Lieb. The next corollary contains 
the main physical and mathematical consequence of the Hohenberg-Kohn theorem 

Corollary 1.4 (HK Theorem and Coupling problem). Let p E 7^(R^) be a measure and ip E 
j^i(R3iV) antisymmetric function with \\ip \\2 = I- There exists a universal map from 

single-particle density p{xi) to a pair densities 72 ( 2 ^ 1 , 3 : 2 ) - or even k-density 'yj^(xi,..., x^) - 
which gives the exact pair density of any N-electron molecular ground state '0(xi,..., x^v) in 
terms of its single-particle density. 


Proof. Consider 

Ip E argmin{(V?, Heiip) : ip ^ p} 

and define by 7 ^ the universal k-point density of that minimizer, i.e. 


7fc(2:i, • 


) = J 


= 


,XN)fdk+i ■ ■■dxN 


□ 
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Of course, 7 ^ may be not unique, since the ground state may be degenerate. The existence 
of a map T : R^ which determines the density of 72 from p is highly not trivial and not 

feasible because still uses high dimension wave function. 


1.4. “Semi-classical limit” and Optimal Transport problem. A natural approach to un¬ 
derstand the behavior of the Hohenberg-Kohn functional F^^[p] ( |rT 2 i i is to study separately, 
the contributions of the kinetic energy and the Coulomb interaction electron-electron energy. A 
well-know method in DFT literature, known as adiabatic connection, is to neglect the kinetic en¬ 
ergy, in which the electron-electron interaction is rescaled with a real parameter A while keeping 
the density p fixed [55l|88] 




i)^p 


y 

T 


N N 

^ + E E 


i=lj=i+l 



(1.15) 


The strictly correlated electron limit (A 00 ), up to relaxing in the space of probability 
measures, was first considered by two papers in physics literature: Seidl m and Seidl, Gori- 
Giorgi &: Savin |83| . In |45| . Gori-Giorgi, Seidl and Vignale interpreted the strong-interaction 
regime as a mass transportation problem with a Coulomb cost. 

Later, an equivalent limit, so-called “semi-classical limit’\ was made mathematically rigorous 
in the two particles case by Cotar, Friesecke Sz Kliippelberg m- They considered the Hohenberg- 
Kohn functional = F^^ as a function of both p and ^ (A = 1) and proved that - up to 

passage of the limit h ^ - the Fhk reduces to the following functional 

F[p] = inf I f - — ^ — :'y 2 {x,y)dxdy:'y 2 £ A 2 ,'y 2 ^ p\, (1-16) 

L Jro \x-y\ J 

where 72 —>■ p means that 72 satisfies the equation and the set A 2 of admissible pair density 

functions is defined by the image of A under the map '0 ^ 72 • As pointed out in [2ZI, instead of 
the corresponding single particle density case, we do not know any characterization of the space 
of admissible pair density function A 2 - 

We state in the next theorem the semi-classical limit for the 2-particles case. The general case 
of N particles is still an open problem. 


Theorem 1.5 (“Semi-classical limit” for N = 2, Cotar-Friesecke-Kliippelberg, [27|)* Let p : R^ —)• 
[0, 00 ] he a probability density such that y/p E H’^(R^), N = 2 and / : R ^ R fee the Coulomb or 
the repulsive harmonic interaction potentials. Then, 

= min ^^ + f{^2 - 

min / f{x 2 - xi)d'y{xi,X 2 ) =: F°^Ip], 
h—^0 7Gn2(IR®,p) Jr3N 

where the set n 2 (R®, p) denotes the set of measures 7 E P(R®) having p as marginals, i.e. (ej){t 7 = 
p where, for 2 = 1 , 2, : R® = R^ x R^ ^ R^ are the projection maps. 

The main difficulty in proving Theorem |1.5| is that for any transport plan given by a density 
0 , 7 = |0(a;i,..., X]s[)\‘^dxi ... x^, it may be that 0 ^ (R^-^), 0 ^ L^(R^-^) and T[0] = 00 . 

Moreover, smoothing the optimal 7 does not work at this level, because this may change the 
marginals of the problem ( |1.17[ |. Cotar, Friesecke Sz, Kliippelberg developed a smoothing technique 
in order to deal with this problem without changing the marginals. A complete proof of the 
previous theorem can be found in m- 

Let us precise somewhat the terminology. In Quantum Mechanics, the “semi-classical limit” 
has a precise meaning: it is an asymptotic regime for the Hamiltonian dynamics of a Quantum 
system defined in a Hilbert space and it is given by a Weyl- Wigner quantization (or quantization 
by deformation studied in the more abstract context of Poisson manifolds). In those specific cases, 
the limit h ^ 0 is called “semi-classical” limit, because the first and second order terms of an 
asymptotic expansion of the Hamiltonian operator is given by “classical” terms, functions of the 
Hamiltonian function in a symplectic manifold ESI EH- 

In DFT context, that limit seems to be up to now merely a question of re-scaling. The 
minimizers of F^^ [p] are candidates of ansatz to develop approximating methods to compute the 
ground state energy of the Electronic Hamiltonian. 

r)p''T' _ OT 

At that point, it is natural to define the DFT-Optimal Transportation ground state Eq 
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^DFT-OT 


[p] = inf |r[V>] + Kie['0] + EotW\ : V’ S 

= inf |rQM[p] + Uie[p] + Eot[p] ■■ P 6 ^|. 


(1.18) 


where p represents the single particle density (see ( |1.5[ |); Tqm are defined, respec¬ 

tively, by 


^QmIp] = inf |t[V>] : V' S "H, i/' ^ p| and 


7 OT 


[P] = 


r N N 

“{LSS 


f{xj - Xi)d'y(xi, ... ,xjv) : 7 S n((M'^)^,p) 


(1.19) 
having p has 


where n((R^)-^,p) denotes the set of probability measures 7 : 
marginals. 

The problem of minimizing (i.igi is called DFT-OT problem. Notice that, in 

the case where / is the Coulomb potential, 

Vee[ 72 ] > Eqt[p], for every probability density p C 

Finally, taking the infimum in both sides of the previous equation, we have 

Theorem 1.6 (Cotar-Friesecke-Kluppelberg, |27p. Consider f{\xj — Xi\) = l/\xj — Xi\. For 
every N, and any potential vr C the density functional with electron-electron 

interaction energy is a rigorous lower bound of the Electronic minimization problem l |1.4| | 


Eq > Eq 


DFT-OT 


It turns out now that the central question of DFT-OT problem with Coulomb potentials is 
to characterize the minimizers 7 of E^'^[p]. Those minimizers could be used as ansatz in order 
to find minimizers of the Hohenberg-Kohn F^^[p]. In [82], the physicist Seidl formulated the 
following conjecture. 

Conjecture 1.7 (Weak Seidl Conjecture [82]). There exists a deterministic minimizer of E^'^[p] 
with Coulomb potential/cost 


N N 

c(3;i,...,3;jv) = 'E E 

i = l 


1 



In other words, for Coulomb-type electron-electron interactions, Seidl conjecture^ state that, 
at least for radially symmetric measures p, there exists an optimal measure 7 for E^^ of form 
7 = (Id X Ti X • • • X Tiv—i)fiP, with Ti,... ,T/v-i : R^- This maps Ti are called co-motion 

functions and, due to the symmetries of the problem, a natural group law is required for them, 
namely that for every i,j there exists k such that Ti o Tj = T^; this is satisfied if for example 
Ti = (cyclical case, see Equation ( |2.9[ |). 

It is natural to ask the same kind of question for other cost functions. In the following we 
mention some interesting results: 


(i) In the 2-electrons case, the minimizers of Eq^^ are well-understood m and there 
exists a map which correlates the density of a given minimization plan of E^'^[p], see 
Theorem nil 

(ii) Gangbo Sz Swi§ch (without being aware of Seidl) reinforces the conjecture 1 1. 7| showing in 

m that for the attractive harmonic cost c(a;i,..., xn) = \^i ~ ^j\‘^ there exists 

a unique optimal plan, and it is a deterministic one. 

(iii) Colombo, De Pascale Sz Di Marino m answer affirmatively (strong) Seidl conjecture in 
the one dimensional case (d = 1 ), as we will describe in section 


^There is also a stronger version of the conjecture, the Strong Seidl Conjecture, where the author descri bes 
explicitly a pos sible optimal map in the radially symmetric case: for the precise statement see the conjecturejS^^ 
and the lemmalS.Slin sectionisl 
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(iv) The authors proved that “Seidl’s conjecture” is n ot tr ue for the repulsive harmonic cost 
c(a:i,. .. ,X]S[) = ~\^i ~ theorem 

for a particular single particle density, there exists 
i.e. an optimal map T which are not differentiable at every point, 
and prove those statements, see section]^ 


4.6 1 . Moreover, for the same cost and 
"Jractal-like” optimal transport maps. 
We are going to state 


2. DFT MEETS Optimal Transportation Theory 

2.1. Couplings and Multi-marginal Optimal Transportation problem. 

Definition 2.1 (Coupling). Let (X, /i.), (T, u) be two probability spaces. A coupling 'y is a measure 
on the product X x Y such that 7 admits p and v as marginals, respectively, on X and Y, i.e. 

(ei)tt7 = M, (e2)t|7 = i' 

where 61,62 are respectively the projection maps (x,y) £XxY\—^xGX and {x,y) £ X x Y 
y£Y. 

It is easy to show that a coupling of p and u is equivalent to 

(i) For all measurable sets A G X and B G Y, one has 'y[A xY] = p[A] and 7 [X x B] = v>[B]. 

(ii) For all integrable measurable functions u, v on X, Y, 

/ {u{x)+ v{y))d‘y{x,y) = / u(x)dp{x) + / v{y)du{y) 

JxxY Jx Jy 

A trivial coupling of (A, p) and (Y, v) is given by the product measure p x u. The symbol 
p ^ tp defined in eq. | |1.9[ l means that the measure 7 = \'tp(xi ,... ,xj^)\‘^dxi .. .dxj^ is a coupling, 
for instance, between pi and pfc-i- In the classical theory of 2 -marginals optimal transport [ 3 | 86 ] 
a coupling is also called a transport plan. 

In the case of a 2-molecular system with Coulomb interactions, up to passage to the limit 
h ^ 0, the theorem |1.5| states that there exists a special coupling 7 , minimizer of ( |1.17[ l, which is 
concentrated on the graph of a measurable function. That motivates the next definition 

Definition 2.2 (Deterministic Coupling). A coupling 7 of two probability spaces {X,p) and 
(Y, u) is said to be deterministic if does there exists a measurable function T : X ^ Y such that 

In that terminology, the main question is to understand when the coupling given by the DFT- 
Optimal Transport l |1.19[ | is deterministic. Moreover, Seidl’s conjecture can be rephrased in the 
following way: Does exist a deterministic coupling in the Optimal Transportation problem 
for Coulomb costs? 

2.2. Multimarginal Optimal Transportation Problem. The “semi-classical limit” in the 
2 particles case (thm. Elf and Seidl’s Conjecture EH motivate us to introduce a more general 
coupling problem. 

We denote by N the number of particles (marginals), d the dimension of the space where 
those particles live (typically d = 3) and I = {1,... ,N} an index set. Take Xi : ^ R as the 

position of a single particle, then the configuration of a system of A’-particles is given by a vector 
X = [xi , X 2 ,..., xtv ) more precisely, by a function x : I ^ R*^. 

Let us consider the distribution of the z-th particle given by a density function pi : R^ [0,1], 
i £ I and c : X R^^^ [0,oo], a continuous cost function. We set pi = piC^,\li £ 7, 

where is the Lebesgue measure on R^, and we want to characterize the optimal coupling of 
the Monge-Kantorovich problem 

(MX) inf / c{xi,... ,XN)d'y{xi,... ,xn), (2.1) 

where , pi ,..., pn) denotes the set of couplings 7 (^ 1 ,..., xjsj) having pi as marginals, for 

all i £ I, i.e. the set n(R^^,/ri,...,/r^v) = {7 G P(R'^^) : {ei)^'y = pi,'^i £ 7} and, for all 
i £ I, ei : R'^^ R^ are the projection on the z-th component and P(R^-^) denotes the set of 
probability measures in R^^. The elements of n(R^-^, pi,..., pn) are called transport plans. In 
the DFT-OT problem it was mentioned that an interesting case is when all marginals are the 
same and equal to a measure p on R^^ (it is quite natural in the DFT framework as the electrons 
are indistinguishable). In that case, we denote the set of transport plans simply by njv(/z). 
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Notice that the existence of a minimum for is not a big deal and it is quite standard 

in Optimal Transport Theory. Indeed, pi,..., trivially non empty, since the in¬ 
dependent coupling pi X • • • X piv belongs to pi,, piv); the set n(R‘^-^, pi,..., piv) is 

convex and compact for the weak*-topology thanks to the imposed marginals; and moreover the 
quantity to be minimized 7 1 -^ f cd'y is linear with respect to p. Hence, we can guarantee the 
existence of a minimum for ( | 2 . 1 [ | by imposing a very weak hypothesis on the cost function, such 
as lower-semicontinuity. 

However, we are interested in characterizing some class of optimal transport plans or, at least, 
understand when that optimal coupling is deterministic. Thus, we would like to study when 
minimizers for (Ad/C) correspond to minimizers the multi-marginal Monge problem 

(Mn) min [ c{x,T2{x), ... ,TNix))dfii{x), (2.2) 

where T/v = {T = ■ Ti is a Borelian map such that T^^fpi = pi,Vi E {1,... ,N} and Ti = 

Id}. As usual, it is easy to see that the set of admissible plans in Sncludes^’ the set of Monge 

transport maps in | |2.2[ |: in fact, given a transport map T = (7d, T 2 ,..., T^v), we can consider 
a measures 77 ^ on (R^^ defined by 77 ^ = (IdxT 2 x ••• x Tiv)tiP- Let h : (R*^)^ ^ R be a 
7 x-nieasurable function, then 

/ h{xi,...,xpf)d'yTixi...,Xff) = h{x-i,T2{xi),... ,Tff{xi))diii{xi)-, 

J (Rd)^ 

and if we have 7 T-nieasurable function / : R^ [ 0 , 00 ] 

Vie I, / f{xi)d'yT{xi,...,XN) = / f{xi)dfj,i{xi), 

jRd 

and so, 7 ^ is a transport plan. In particular, the value of the minimization over Monge-Kantorovich 
couplings is less than or equal to the value of the minimization over the deterministic couplings 


min(Ad/C) < inf(Adiv)- 

Contrary to in the existence of the optimal deterministic coupling in is not obvious and 

the difficulties do not lie on the “multi-marginality” of the problem. In fact, at that point, it is 
not possible to apply standard methods in Calculus of Variations even in the classical 2-marginals 
case in order to guarantee the existence of the minimum \2.2 \ . 

Intuitively, the difference between ( |2.1[ | and \2.2\ is that in almost every point x\ E R*^ is 
coupled with exactly one point xi E R^ for each i E {2,..., V} whereas throughout the problem 
we allow the splitting of mass in the “transportation” process. 

The problem ( 23 . when N = 2, was studied by L. Kantorovich in 1940’s 53 Eg , without being 
aware of the famous article of Monge m- In 1975, L. Kantorovich won the Nobel Prize together 
with C. Koopmans ‘for their contributions to the theory of optimum allocation of resources”. 

An important disadvantage in using the relaxed approach in the multi-marginal OT is 

that the set of optimal transportation plans could be very large and we could need to select some 
special classes of transportation plans. 

In the 2-marginal setting, the two formulations and \2.2\ are proved to have the same 

value for general Polish Spaces by Pratelli in ED- 


Theorem 2.3 (Monge equals Monge-Kantorovich in the 2-marginals case, ED). Let{X,p),{Y,u) 
be Polish spaces, where p E P{X) and u E P(V) are non-atomic probability measures. If c : 
A X y —)• R U { 00 } is a continuous cost, then 

min < / c{x,y)d'y{x,y) : ■y £ Il{X x Y, p, u) 

[ JxxY 

= inf{^c(x,T(.;))dMx), : f f y Borel } 

where n(A x y) denotes the set of 2-marginals transport plans ■y G X x Y having p and u as 
marginals. 

In the particular case of the multi-marginal OT problem in and 113 when all marginals 

are equal to pC^, we can apply the previous theorem ( |2.3[ | on the Polish spaces X = R^,y = 
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1 measures = pC^ and u = (e 2 , •.., en)ji 7 , for every 7 E n(R^^, pC^), and obtain 

the following corollary 


Corollary 2.4. Let pi = ■ ■ ■ = = pC^ be probability measures on with density p and 

c : (R^)^ ^ R fee a continuous cost function. Then, 


min(A4/C) = inf(AlAr). 

Notice that , in general, the equivalence between and \2.2\ is not an immediate consequence 


of the theorem 


2.3 


In particular, remark that the image measure of y = is not prescribed, 

but only its marginals. 

Finally, we conclude this section by giving an example of a particular continuous cost function 
where the optimal coupling 7 is not deterministic, even in the 2 particles case. The nonexistence 
of an optimal transport maps typically happens when minimizing sequences of 7 exhibit strong 
oscillations. 


Example 2.5 (G. Carlier). Suppose p uniformly distributed in [0,1], u uniformly distributed in 
[—1,1] and the continuous cost c : R^ ^ R given by 

c{xi,X2) = {xl - 

Since c > 0 we have that an optimal plan for this cost is given by 7 = + since 

j cd 7 = 0. We want to show there doesn’t exit a map T such that the cost attains zero: suppose 
in fact that j{x^ — T{x)‘^)‘^ dx = 0 but then one should have T{x) E {a:, —x} almost everywhere. 
Hence, 

T{x) = {XA - X[o,i]\a)^ 

for a measurable A C [0,1]. But then Tjj// = C^\a + ^^\a ^ where a: E A iff —x E [0,1] \ A, 
and so we obtain a contradiction since T wouldn’t admissible. It is however easy to construct 
admissible maps Tn, linear on [i/2n, (i + l)/ 2 n), such that f c(xi, Tn(xi))dp(xi) tends to zero. 

Remark that, we state by simplicity the example in the 2-marginal case, but the same conclu¬ 
sion and arguments apply in the larger interactions case under the cost 

N N 

c(xi,...,xn) = '^ ^ 

i=lj=i+l 



Figure 1 . (example \2.5\ ) Left: support of the optimal coupling 7 . Right: 
optimal coupling 7 . The simulation has been performed by using the iterative 
method described in section[El 


2.3. Dual formulation. We consider the formal convex dual formulation of 113 


{Kn) sup 


N N 

E Ui{xi)dfii{xi) : E Ui{Xi) < c{xi, 

i = l i=l 



(2.4) 


where Ui E L^(R^,pi), called Kantorovich potentials or simply potentials, are upper semicontin- 
uous for all i E {1,..., W}. This is often called the dual problem of 1 ^ . 
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Notice that, when all marginals are equal to fii = p{x)C^ and the cost is symmetric, we can 
assume the Kantorovich potentials Ui{xi) are all the same u{xi), and we can rewrite the constraint 
in | |2.4[ | as 

N 

E u{Xi) < c(xi, ... ,xn)- 

i = l 

Among all Kantorovich potentials, a particular class is of special interest. 

Definition 2.6 (c-conjugate function). Let c : (R^)^ Ru{oo} a Borel function. We say that 
the N-tuple of functions {ui,... is a c-conjugate function for c, if 

Ui{xi) = inf < 0 ( 3 : 1 , • • ■ ^ ^ j ^ f 5 yi = 1,..., N. 


As we will see in the next theorem |2.7| the optimal potentials ui,...,Ujv in {K,n) c- 
conjugates, and therefore semi-concave. So, they admit super-differentials dui{xi) at each point 
Xi G R^ and then, at least in the compact case, we can expect that for each Xi G R^, there exists 
(xj)j^i G such that 

Ui{xi) = c{xi,. ■. ,xn) - y^^Uj{xj). 

Moreover, if c : (R^)^ ^ R is differentiable and \ui{xi)\ < 00 , for some Xi G R*^, Ui is locally 
Lipschitz. The well-known general result linking both IMJ and l |2.4[ | has been proven by Kellerer 
in 1984. 


Theorem 2.7 (Kellerer, |48p. Let (Aij/ii),..., (A’Ar,/ijv) Polish spaces equipped with Borel 
probability measures /ii,... , /ijv • Consider c : Xi ... Ajv ^ R a Borel cost function and assume 
that c = supxc < 00 and c = infx c > — 00 . Then, 

(i) There exists a solution 7 of the Monge-Kantorovich problem iHl and a c-conjugate 
solution {ui,U 2 , ...,un) of the dual problem (2-4]) ■ 

(ii) “Duality holds”, 


inf 


cd'y : 7 G n(/ii,...,/^iv) ^ = sup 


J Xi...Xj^ 

where D is the set of functions ui : ] 


r ^ r 

■iu 


Ui{xi)d^ii ■.Ui&'D), 


(2.5) 


i = 1 ,..., N, such that 
N 


'^Ui(xi) < c(xi,... ,xn), and -^c - (c - c) < "^Ui{xi) < ^c. 

i=l i=l 

(hi) For any solution 7 of 113, any conjugate solution of \2A\ and 7 a.e. ( 3 : 1 , ..., xx), we 

N 

'^Ui{Xi) = c(xi,. . .,X]sj). 


have 


We remark that despite its generality, Kellerer’s theorem can not be applied directly to 
Coulomb-type costs, since Coulomb cost does not satisfies the boundedness hypothesis c = 
supxc < 00 and can not be bounded by functions in L^(R^, pC^). 

A generalization of Kellerer’s theorem due to M. Beiglboeck, C. Leonard and W. Schachermayer 
[ 6 ], extends this duality for more general costs, but also cannot be applied in this context. 

Recently, De Pascale 1301 extended the proof of duality which can be applied also for Coulomb- 


type costs, see Theorem 3.1 In fact, suppose p = p{x)dx a probability measure on 
not give mass to sets of cardinality less than or equal to d — 1 , then 

N N ^ 

E E 

1=1j=i+l 


mm 

7enjv(M) 


/ 


-dy = max 

II 


[ 


Nu{x)p{x)dx 


which does 


(2.6) 


where D is the set of potentials u G L^(R^, p) such that 


N 

E“' 

i=l 


(a:*) < 






iSifLiP — almost everywhere. 


(2.7) 
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Moreover, we get the conclusions (z) and {iii) of the theorem |2.7| Physically, the constraint 
means that, at optimality, the allowed Nd configuration space should be a minimum of the 
classical potential energy. 

Remark 2.8. More generally, the main theorem in EQl shows that the Kantorovich duality holds, 
for instance, for costs of form 


N N 


E E 

i=l j = 



S> 1. 


2.4. Geometry of the Optimal Transport sets. Our goal in this section is to present the 
relation between the support of a coupling 7 and optimality in the Monge-Kantorovich problem 
In the following, we will just summarize key results necessary to present recent developments 
of optimal transportation with Coulomb and repulsive harmonic-type costs. 

Roughly speaking, it is well-known in optimal transport theory with 2-marginals that, for a 
wide class of costs c, a coupling 7 is optimal if, and only if, the support of 7 is concentrated in a 
c-cyclically monotone set ElCElllO]: 

Definition 2.9 (c-cyclically monotone, 2-marginal case). Let X, Y be Polish spaces, c : X xY 
R U { 00 } be a cost function and T be subset of the product space X x Y. We say that F is a 
c-cyclically monotone set if, for any finite couples of points {{xi^yi)}^^ C P, 

M M 

E c{xi,yi) < E 

i=l i=l 


In |67| . Pass suggested a possible extension of the concept of c-cyclically monotone sets for the 
multi-marginal case. For our proposes, it will be enough to consider this notion in 

Definition 2.10 (c-cyclically monotone set). Let c : > R be a cost function. A subset 

r C R"^^ is said to be c-monotone with respect to a partition {p,p‘^),p G I, if, for all x = 
{xi,.. .,Xff),y = {yi,.. ■,yN) e r 

c(x) -I- c{y) < c{X{x,y,p)) + c(Y{x,y,p)), 

where X{x,y,p),Y{x,y,p) £ R^ are functions obtained from x and y by exchanging their coor¬ 
dinates on the complement of p, namely 


Xi{x,y,p) = 


Xi, if i G p 
yi, ifiGp'^ 


and Yi(x,y,p) = 


yi, ifi&P 

Xi, ifiGp'^ 


MiG I 


We say that F C R'*^ is c-cyclically monotone, if it is c-monotone with respect to any partition 

pel. 


The notion of c-cyclical monotonicity in the multi marginal setting was studied by Ghoussoub, 
Moameni, Maurey, Pass |42l 14311671169| and successfully applied in the study of decoupling PDE 
systems [44| . For Coulomb costs, it was useful in order to characterize the optimal transport maps 
in the one-dimensional case m- 

More recently, Beiglboeck and Griessler [Jj, presented a new notion called finistic optimality 
inspired by the martingale optimal transport, which is equivalent, in the 2 -marginals case, to the 
definition im 

The next proposition gives a necessary condition on the support for an optimal transport plan 
7 . For the proof, which is based on the 2-marginal result m, we refer the reader to Lemma 3.1.2 
of [67]. 

Proposition 2.11 (Support of transport plan are c-cyclically monotonic). Let c : (R^)-^ [0, 00 ] 

be a continuous cost and , p]s[ absolutely continuous probability measure on M.^. Suppose 

that 7 E n((R^)-^, ,. .., pjv) is an optimal transport plan for Multimarginal Monge-Kantorovich 

Problem and assume (A4/C) < 00 . Then 7 is concentrated in a c-cyclically monotone set. 
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2.5. Symmetric Case. We are going to show in this section that we can reduce our study 
of the Monge-Kantorovich and the Monge problems ( |2.2[ ), respectively, to a suitable class of 

transport maps (called cyclic maps) and transport plans (symmetric plans). Having in mind the 
DFT-OT problem, it is natural to suppose that the cost functions does not change if there are 
some symmetries among the configuration system of particles or, equivalently, we can suppose the 
couplings 7 G TlN(fi) of those particles are cyclically symmetric. 

Definition 2.12 (Cyclically symmetric probability measures). A probability measure 7 G 
is cyclically symmetric if 

[ f{xi,...,XN)d'y=f /(ct(xi, ... ,a;iv))<i7> V/e C(R^‘*) 

jRdN J-RdN 

where u is the cyclical permutation 

a{xi,X2, . . . ,Xn) = {X2,X3, • ■ • 

We will denote by nj^^(/ii), the space of all symmetric probability measures on having the 
same first marginal pi. 


Notice that if {pi,..., pjsj) are the marginals of 7 G , pi), then all marginals 

pi,..., PN are necessarily the same: V / G C(R‘^), 


Wiel, / f{xi)dpi{xi) 
JR^ 


where f{xi, ... ,xjv) = fi^i)- Briefly, 

(ei)tt7 = G I, 

where is the canonical projection. 

As previously, we can define the multi marginal symmetric Monge-Kantorovich Problem as 
follows: 

{MlCsym) inf / c{xi,... ,XN)d-f{xi,. .. ,xn), (2.8) 

and the associated multi marginal cyclically Monge Problem 

(Mcycl) inf f c(a;,r(i)(a:),r(2)(x),...,r('^-l)(a:)), (2.9) 

= J]^d 

jpiN) _j 

Notice that, (Aicyci) is related to (Ad) by choosing T^p = p, Ti+i(x) = T^'^^{x) for every i G 
- l},rW = Id. 

Proposition 2.13. Suppose p S is an absolutely continuous probability measure with 

respect to the Lebesgue measure and c : (W^)^ M. is a continuous symmetric cost. Then, 

■m{(MK) = miiMlCsyrr,). 


'■ f f{a''-~’-(xi,... ,XN))d'y 

: / f{x-l)dpi{xi), 


Proof. The fact that the infimum of is less than or equal to | |2.8[ l is obvious. Now, we need 
to show that for every transport plan 7 G we can associate a symmetric plan in n^’^(/i). 

Indeed, given 7 G TIn{p), we define 



'^# 7 - 


and finally we notice that 7 has the same cost as 7 . 


□ 


We remark that prove m.i{fAcycl) ^ inf(Adiv) is not obvious . However, we can avoid this 
problem in order to prove the equivalence between the Multi-marginal cyclic problem and Multi¬ 
marginal OT problem, by noticing that inf(Adsym < inf(AdcycO- "^^e last part was proved by M. 
Colombo and S. Di Marino m- We give the precise statement of this theorem. 
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Theorem 2.14 (M. Colombo and S. Di Marino, [24p. Let c : (R^)^ M. be a continuous 
function and bounded from below. If /i has no atoms, then 

{■AdtCsym') — (-^cycz)* 


Finally, we have the equality between the Multi-marginal cyclic problem and Multi-marginal 
OT problem. 

Corollary 2.15. If // has no atoms and c : (R^)-^ is a continuous function and bounded 

from below, then 

inI{MlCsyin) = mi{MK.) = inf(AlAr) = inf(Alcycz)* 

The last result of this section is due Pass [70], which states that optimal transport plans for 
symmetric costs can be supported on sets of Hausdorff dimension equal or bigger than 2N — 2. 

Theorem 2.16 (Pass, |7Qp. Let c : (R*^)-^ ^ R 6 e a symmetric cost and radially 

symmetric and absolutely continuous with respect to Lebesgue measure in R^. Suppose that for 
every radii (ri, ...r^) the minimizers 

{xi,, xn) E argmin|y^c(yi, y2, • ■ ■ Pn) 'f^ot all co-linear. 

Then, there exists solutions 7 for the Monge-Kantorovich symmetric problem l |2.8[ l whose support 
is at least (2N — 2)-dimensional. 


3. Multi-marginal OT with Coulomb Cost 

This section is devoted to the summary of the results present in literature on the multi 
marginal optimal transport with Coulomb cost. We already highlighted the fact that the inter¬ 
esting case is when all the marginals are equal, since this is the physical case. We point out that, 
as always, the focus will be on characterizing the optimal plans as well as looking for cyclical 
transport maps: in this direction an important conjecture has been made by Seidl in its seminal 
paper [82], where he gave an explicit construction of a map for measures with radial symmetry 
in R^. In the last years the conjecture was proven to be right: first in the ID case when p is 
the uniform measure on an interval, and then generalized to any diffuse measure on the real line. 
However quite recently various authors disproved the conjecture in the case R*^ for d > 2 and for 
TV = 3 marginals, yet the Seidl map is still optimal for some class of measures. 

We begin by analyzing the general problem, and then we will proceed to look at the radial 
case. 


3.1. General theory: duality, equivalent formulations and many particles limit. We 

begin by stating the duality theorem in this case proved by De Pascale [301, by means of F- 
convergence of finite dimensional problems for which the duality is classical. 

Theorem 3.1. Let p E P(R'^) be a measure that is not concentrated on a set of cardinality 
smaller or equal than TV — 1. Then the duality ( |2.6| | holds and the dual problem has a bounded 
maximizer. 


The boundedness of the maximizer u can let us prove also some regularity properties. The 
first one is the fact that any optimal plan is supported away from the diagonals, while the second 
one proves second order regularity of u. 

Lemma 3.2 (Corollary 3.13 in [30]). Let p as in Theorem \3.1\ then there exists o > 0 such that 
for every optimal plan 7 we have supp( 7 ) C {{xi, ..., xn) : \xi — Xj \ > a, \/ i ^ j}. 

Lemma 3.3. Let p as in Theorem \3.1\ and u be a bounded maximizer in the dual problem in 
| |2.6[ l,' then u has a p-representative that is Lipschitz and semi concave. 


Proof. As in the classical case we first prove a structural property of the maximizer u, namely 
that the constraint 113 is saturated, that is (u,... ,u) is a c-conjugate TV-tuple. Let us consider 


w(xi) = 


ess-inf 


l<i<j<N I * i=2 


where the ess-inf is made with respect to the measure p'^^ It is obvious that w{xi) > u(a:i) for 
p-SL.e xi, thanks to { 23 . Suppose that ic > u in a set of positive measure; but then {w, u,.. . ,u) 
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would be a better competitor in the (not symmetric) dual problem contradicting the fact that 
u is a maximizer (we use the fact that the symmetric dual problem and the not symmetric one 
have the same values). Now we have that w is defined everywhere and so we can talk about its 
regularity. 

Let d > 0 be a number such that there exists points pi,...pN in the support of p such 
that \pi — Pj\ > d; this number exists as long as ttsupp(/i) > N. Let us fix e > 0 such that 
^ ^ + 2 N||u||oo and e < d/ 8 ; in this way it is true that if we define 


= ess-inf 

|a:i—a;o|>£, i=2,...,iV 


E 

l<i<j<N 



N 1 

i=2 ) 


we have Wxq = w on B(xo,£). In fact in the definition of w we can choose xi among the pi (or 

very close to them, as they belong to the support) such that \xi — Xj\ > d/2 — 2e > d/4 for every 

i 7 ^ j, and so we have that w < _|_ NjlnUoo but then it is clear that for any ( 0 : 2 , • ■ •, xj\f) 

such that \xi — xq\ < e we will have 


y; -— -- - y > — > lo(xi) VxieB(xo,e); 

l<i<i<N I*' 

this proves that in fact Wxq = to on the set B^xq^s) and so in particular also in B(xq^£/2). But 
in this set Wxq is Lipschitz and semi concave since it is an infimum of uniformly functions on 
B{xq,£/2). Moreover the bounds on the first and second derivatives don’t depend on xq but only 
on £•, that is fixed a priori, and so by a covering argument we obtain the thesis. □ 


Another interesting reformulation of the Coulomb-like problem, or more generally when we 
have only interaction between two particles, can be found in m where, seeking a dimensional 
reduction of the problem, the authors use the fact that if 7 E njv(/^) is a symmetric plan then 

f E I- - - I ^ . d(ei,e2)i'r, 

l<i^<N ^ 2 / jR 2 d \x-y\ 

that is, since the electron are indistinguishable, it is sufficient to look at the potential energy of 
a couple of electrons and then multiply it by the number of couples of electron. It is clear that 
72 = (^ij ^ 2 )tt 7 is a 2-plan whose marginal are still p\ we will say that rj E Tl 2 {p) is A’-representable 
whenever it exists 7 E niv(M) such that rj = (ei,e 2 )tt 7 . The equivalent formulation they give is 


mm 


I 

J-gNd 


E 


\xi - Xj\ 


dy = 


mm 

T7Gn2(Ai) 

7 } is A^-representable 


nf 


2 ) J^ 2 d \x - y\ 


dp. 


Unfortunately the conditions for being 7V-representable are not explicit and they are very 
difficult to find; this is correlated to the A’-representability for matrices, but here, since we are in 
the semiclassical limit, we deal with densities instead. In m the authors propose, as a method 
for reducing the dimension of the problem, to substitute the condition of being A-representable 
with that of being /c-representable with k < N; the resulting problem will give a lower bound to 
the SCE functional . 

We present also a general theorem regarding the many particles limit, that embodies the well- 
known fact that when N ^ 00 then the solution is the mean field one 


Theorem 3.4 f|28p. Let p E P(R'^). Then we have that 

lim min [ -- - -- d'y = [ -—-—- dp (g) p. 

( 2 ) ^^^N(n)jRNd \xi -Xj\ J^2d \x-y\ 

In terms of DFT we are saying that E^'^{p) = (^)_E^^(p)-|-o(A^), where is the normal¬ 
ized mean field energy E^^{p) = f dp(x) dp{y). 

Remark 3.5. The statement about the Coulomb cost in the physical case is quite classical. In 
fact, for a measure p E L'^'^^(R^) the Lieb-Oxford bound holds [561158| : 

(^']e^^{p) > E°'^(p) > N^E^^^ip) - [ p‘^/^{x) dx, 

V 2 / 

and so the conclusion is immediate. However, in m the proof is completely different and relies 
on the fact that a measure 7 E Tl 2 {p) that is A-representable for every N must be in the convex 
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envelope of measures of the type i/ (g) z/, and then on a direct computation using the Fourier 
transform. In particular, aside from the Coulomb cost in dimension 3 they prove the theorem 
for more general costs with binary interaction of the form ~ ^j)-> where I G 

n L^(R^) satisfies l{z) = l{—z) and also 1 > 0 . 

3.2. The Monge problem: deterministic examples and counterexamples. In this section 
we will illustrate the cases in which we know that there exists a deterministic solution; in those 
cases moreover there is also a result of uniqueness. We remark that these are the two extreme 
case, namely the case N = 2 (and any d) and the case d = 1 (and any N). We will just sketch the 
proofs in order to make clear the method used here, and why it is difficult to generalize it. We 
begin with the 2 -particles case, in every dimension d: this result was proved in by means of 
standard optimal transport techniques. 

Theorem 3.6. Let /r G be a probability measure that is absolutely continuous with respect 

to the Lebesgue measure. Then there exists a unique optimal plan 'yp G 112 (/i) for the problem 

min f - -- c/ 7 . 

7en2(Ai)7K2d \x\ — X 2 \ 

Moreover this plan is induced by an optimal map T, that is, 7 = {Id, T)^p, and T{x) = x+ 
p-almost everywhere, where (p is a Lipschitz maximizer for the dual problem. 


Proof. Let us consider 7 a minimizer for the primal problem and 0 a bounded and Lipschitz 
maximizer of the dual problem (it exists thanks to Lemma |3.3| l. Then we know that 

F{xi ,X 2 ) = 0 ( 3 ^ 1 ) + 0 (^ 2 )-^ 0 for p (g) /d-almost every {x\,X 2 ), 

FI - ^ 2 \ 

and we know also that F’ = 0 7 -almost everywhere. But then F has a maximum on the support 
of 7 and so VF = 0 in this set; in particular we have that V0(a:i) = — , on the support of 

\x\ —X2 u 

7 . Finding X 2 we have that X 2 = xi — (^i) “ ^(^ 1 ) the support of 7 and this implies 

7 = (7c/, as we wanted to show. □ 

The first positive ^/-marginal result for the Coulomb cost is instead shown in m where, in 
dimension d = 1, the authors can prove that for non-atomic measure an optimal plan is always 
“induced” by a cyclical optimal map T. 


Theorem 3.7. Let p G F(M) be a diffuse probability measure. Then there exists a unique optimal 
symmetric plan 70 G Tl^^^ip) that solves 


mm 


/ 


E 

l<i<j<N 


1 

\xj - Xi\ 


c/ 7 . 


Moreover this plan is induced by an optimal cyclical map T, that is, 7 q = 
where 77 ^ = (7c/, T, , • • •, T^^~^'^)i^p. An explicit optimal cyclical map is 

^ ('F,7^(F,.(x) + l/iV) ifF^ix) < iN-l)/N 

yF^^{Ffj,{x) 1 — 1/N) otherwise. 

Here F^{x) = p{—oo,x] is the distribution function of p, and F/T^ is its lower semicontinuous 
left inverse. 


Proof. We begin by observing that if 7 is a symmetric optimal plan then a stronger statement 
than Proposition |2. 1 1] holds, namely we have that for every x,y ^ supp( 7 ): 

c{x) -h c(y) = mm{c{X{x,a{y),p)) -h c{Y{x,a{y),p)) : Vp C {1,... ,/V},Vcr G 

where c is the Coulomb cost on the 77-tuple of points x = {xi,... ,xn) and cr acts on the indices. 
The key point here is that one can show (Proposition 2.4 ED) that this property holds if and 
only if a more geometric condition holds true for x and y: they are well-ordered. This property 
amounts to the fact that, up to permute the coordinates of both points, we have 

XI < yi < X2 < ■ ■ ■ < XN < VN, 

or the other way around. Once this property is proven the rest of proof is rather straightforward, 
since we proved that the support of (ei, e 2 )tt 7 |D? where D = {xi < Xj and yi < yj, Vz < j}, is 
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monotone in the usual sense. We remark that the key property does not easily generalize to d > 1 
since it (heavily) uses the ordered structure of R. □ 


Let us notice that in this case we don’t have uniqueness of optimal plan neither of cyclical 
optimal maps, as pointed out in Remark 1.2 in m- 

For the general case we thus have to assume d > 1 and N > 2\ from now on we will look at 
the case in which the measure has radial symmetry. We can easily see that in this case we can 
reduce our problem to a multi marginal problem in R. 


Lemma 3.8 (Radial case). Let E 7^(R^) be a radially symmetric measure. In particular p is 
determined by pr = \ • \^P- Then for every optimal plan 7 E njv(M)j if we consider ■yr = Tt 0 , 
where R : (xi,... ,a:jv) (l^il? • • • ? I^iv|)j we have that yr is an optimal plan for the ID multi 
marginal problem 

min / ci(ri,... ,riv) dr?. (3.1) 

T/en^vCMr) Jrn 

where ci is the reduced Coulomb cost 


ci(ri,... ,rjv) = min 


E 



\xi\ 



Moreover Ei<i<j<jv = /r^ ci(ri,... rjv) d'Tr- 


Given the one dimensional character of this problem, one could expect that the solution is 
similar to the one depicted in Theorem |3.7| In fact in m the authors conjecture a similar 
structure: 


Conjecture 3.9 (Strong Seidl Conjencture). Let p E 7^(R'^) be an absolutely continuous measure 
with respect to the Lebesgue measure, with radial symmetry, and let pr = \ ■ Let 0 = ro < 
ri < ... < ‘Tjv —1 < f'N = such that the intervals Ai = [ri,ri+i) have all the same radial 
measure pr{Ai) = 1/A^. Then let F{r) = pr{0,r] be the cumulative radial function and let 
S : [0, 00 ) [ 0 , cxo) be defined piecewise such that the interval Ai is sent in the interval in 

an anti monotone way: 


S(r) = F-\2i/N - F(r)) ifn^i < r < n 
I F~^(F(r) + 1/N — 1) if H is even 
[F-^{1- F{r)) ifN is odd, 

Then S is an optimal cyclical map for the problem 


F(Sir)) = 


and i < N 
ifrN-i < r < Tat. 


However in the recent papers [251132| this conjecture is proven to be wrong, looking at radial 
measures p concentrated in some thin annulus {1 — 5 < |a^| < 1 + <5}; in particular in |32| the ID 
problem is proven to be equivalent, when 5 ^ 0 , to the repulsive harmonic one, for which high 
non-uniqueness holds and the Seidl map is not optimal (see the following section). However in 
m also a positive example is found, namely a class of measures for which the conjecture holds. 

These results show that the solution to the multi marginal problem with Coulomb cost is far 
from being understood: in particular there is no clear condition on the marginal p for the strong 
Seidl conjecture to hold. 


4. Multi-marginal OT with repulsive Harmonic Cost 


This section is devoted to the study of the repulsive harmonic cost. In DFT-OT problem (see 
section Q, we replace the electron-electron Coulomb interaction by a weak force which decreases 
with the square distance of the particles. The term weak comes from the fact that, in the repulsive 
harmonic cost the, interaction cost function has value zero when the particles overlap, instead of 
00 in the Coulomb case. 

More precisely, we are interested in characterizing the minimizers of the following problem 




mm 


/ 

.'fR 


.Mjv) 


-\xj - Xifd-yixi,... ,xn), (4.1) 

Ui=l 


where pi,, p^ are absolutely continuous probability measures in R^. 

From a mathematical viewpoint this cost has some advantages compared to the Coulomb one, 
since here we can do explicit examples. We can interpret the solutions of l |4.1| l as an ansatz for 
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DFT problem (Hohenberg-Kohn functional ( |1.12[ l) for particles interacting under the repulsive 
harmonic potential. From a technical aspect, this could be an interesting toy model to approach 
the case of Coulomb cost. From applications, some minimizers of l|4.1|l se em to have no particular 
relevance in physics because, as we will see in the examples |4.10| |4.TT] and |4.12| below, certain 
optimal 7 of | |4.1[ | allow particles to overlap. 

We will see that this problem has a very rich structure, that is very different from the classical 
2-marginal case. First of all we notice that minimizers of this problem are also minimizers of the 
problem with the cost c(xi,... ,Xn) = \xi . -\-Xn\^', in fact we have that 



\xf dlMi, 


but this last additive term depends only on the marginals and not on the specific plan 7 . The 
cost c is very particular since it has a wide class of “trivial” optimal plans, that is the ones that 
are concentrated on -h ... + =0; however the structure is very rich, see Lemma |4.3| 

Concerning the existence of minimizers of ( |4.1[ |) we will assume that the measures have 
finite second moments; then existence follows immediately from the equality 


^ N N 

argmin / - Iv 

i=lj=i+l 


\^d'y = argmin 
7 


/ixi 


■ xjv|^rf7- 


In particular, Corollary |2.15| holds in this case too. Notice that the fact that the repulsive harmonic 
cost is smooth and has linear gradient does not make the Multi-marginal Optimal Transportation 
problem easier compared to the Coulomb cost. In fact, in this case the problem is that if we write 
down the optimality conditions for the potentials, in the case described in Lemma |4.3| we simply 
find the condition ici -|- ... -h x^ = const, since in this case the potentials are linear functions. 

However, we can enjoy the symmetries of this problems and build easy Monge solutions for 
some particular cases of ( |4.1[ l, see examples |4.9[[4.10[|4.11||4.12| and EH below. 

Before stating the main result in the multi-marginal setting, we start analyzing the problem 
| |4.1[ | in the 2-marginals case, where everything seems to work fine, just as in the square distance 
case. 


Proposition 4.1. Let fi,u £ P(R'^) and Cuj{x,y) = —\x — y\^ be the opposite of the square- 
distance in R^. Suppose that p is an absolutely continuous with respect to the Lebesgue measure 
in R^ and u has no atoms. Then, there exists a unique optimal transport map T : R*^ R'^ for 

the problem 


min 


I j Ou,{x,y)d') 


7 e n2(iR2‘' 



Cw{x, T{x))dp 



Moreover, T = Vcf), where <j) : R^^ —^M.isa concave function and there exists a unique optimal 
transport plan 7 , that fs 7 = (IcZ x T)^p . 


Proof. This result is a easy consequence of the Brenier’s theorem. Indeed, it is enough to verify 
that, taking C = 2 f \x\'^ dp -\-2 j \x\'^ du, we have 

C-h inf [ —\x — T{x)\^dp{x) = inf [ \x — G{x)\^dp{x) 

G^tJ.=u J 

where G = —T and u = [—Id)'^v. Then by Brenier’s theorem there exists an unique optimal map 
G which can be written as G{x) = ^/ijjix), and 1 /) : R^^ ^ R is a convex map. In other words, T is 
a gradient of a concave function. □ 


Notice that if we suppose p = f{x)dx and u = g{x)dx are probability measures with densities 
f,g concentrated, respectively, in convex sets and assume there exits a constant A > 0 

such that A < f,g < 1/A. Then, T is a (7^’“ function inside 0.f [131114) . 

Example 4.2 (2 marginals case, uniform measure in the d-dimensional cube). Suppose that p = 
1 / = In this case, we can verify easily that the optimal map T : [0, l]*^ —>■ [0, l]*^ is given 

by the anti-monotone map T{x) = (1,..., 1) — x. 

Surprisingly, the next theorem says that we can not always expect Caffarelli’s regularity for the 
optimal transport maps for the repulsive harmonic cost with finitely many marginals N > 2 even 
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when the support of the measures has convex interior (see corollary |4.8| below). Before stating the 
theorem we will prove a lemma that characterizes many optimal plans: 

Lemma 4.3. Let be probability measures on and h : R^ 'R be a strictly convex 

function and suppose c : ([0, 1 ]^)^ ^ R a cost function of the form c(a:i,. .. ,xn) = h{xi + 

• • • + xn). Then, if there exists a plan 7 G ..., concentrated on some hyperplane of 

the form xi . xn = k, this plan is optimal for the multi marginal problem with cost c and 
supp( 7 ) C {x\ + ... + = k} is a necessary and sufficient condition for 7 to be optimal. In 

this case we will say that 'y is a flat optimal plan and {pi}fLi is a flat N-tuple of measures. 


Proof. First of all we show that k is fixed by the marginals pi, and it is in fact the sum of 
the barycenters of these measures. Let Ci = f xdpi\ then let us suppose that there exists 7 G 
T{pi,..., pj\i) that is concentrated on {xi + ... + = k}. Then, using that the 2 -th marginal 

of 7 is Pi, we can compute 

„ N „ N „ N 

k = (xi xn) ^7 = XI / d7 = XI / xdpi = '^Ci. 

i=l i=l i=l 


In particular we notice also that for every admissible plan 7 we have f {xi xn) d'y = k and 

so by Jensen inequality we have 



-h ... + xjsi) d7 > 


h 


(^J (3:1 -f ... + x^) c/7^ = h{k) = j h{xi -h . .. + x^) d'y. 


This proves that 7 is an optimal plan. Thanks to the strict convexity of h, this shows also that if 
7 is optimal then 7 -a.e. we should have xi + x^ = k. □ 


This reveals a very large class of minimizers in some cases, as we will see later. However not 
every A^-tuple of measures is flat: it is clear that we can have marginals such that there is no plan 
with such property: 

Remark 4.4. Let N = 3 and pi = p for every i = 1, 2,3 with p = 1^1 + 1^2 with ui = (—Id)i^i '2 
and ui concentrated on [2, 3]. Now it is clear that the barycenter of /i is 0 but for every 3 points 
xi, X 2 : xz in the support of p we cannot have xi -\- X 2 -\- xz = 0 : two of them have the same sign, 
let’s say xi and X 2 , but then we have \xi + X 2 \ ^ 4 > 3 > which contradicts xi -j-X 2 = —X 3 . 
So there is no an admissible plan concentrated on Hq = {xi -\- X 2 -\- xz = 0}; in fact we showed 
that for every 7 G n(/ 2 ) we have supp( 7 ) D Hq = 0. 

Remark 4.5. In the case N = 2 the condition for which there exists an admissible plan 7 concen¬ 
trated on an hyperplane of the form xi-\-X 2 = k implies that the two measures p and i> = {—Id)^^ 
are equal up to translation. This condition is in fact very restrictive. However when we think at 
the DFT problem and in particular we assume that p = u then the said condition amounts to 
have that p is centrally symmetric about its barycenter and in the context of density of 2 electrons 
around a nucleus this seems a fairly natural condition. 

Theorem 4.6. Let pi = p = ^^^,V 2 = 1,...,A^ be the uniform measure on the d- 

dimensional cube [0,1]^ C R^^, h : R^ R a convex function and suppose c : ([0,1]^)'^ R 
a cost function such that c{xi ,... ,Xjv) = h{xi -!-••• + ^n)- Then, there exits a transport map 
T : [0,1]^ —>■ [0,1]*^ such that T^(x) = x and 

min / cd'y = min j c(x,T(x),... ,T^^~^\x))dp 

Moreover, T is not differentiable at any point and it is a fractal map, meaning that it is the 
unique fixed point of a “self-similar” linear transformation T acting on L°°([0,1]^; [0, 

Proof. We express every 2 : G [0,1] by its base-A^ system, 2 = with a^, G {0, \ ... ,N — 1}. 

Consider the map given by S{z) = ^ ? where S is the permutation of N symbols such 

that tS( 2 ) = 2 + 1, V {2 = 0,1,..., N — 2} and S{N — 1) = 0. A straightforward computation shows 
that 

^ N 

i=l 


(4.2) 
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Let T : —)■ be a map defined by T{x) = T(zi,..., Zd) = {S{zi),... , S{zd)) and denote 

by (a?) = (-^ 1)5 • • • ? i^d))i j = 1,..., N — 1. We will first show that S is a measure¬ 
preserving map. In fact, we can show that there exist functions Sk : [0,1] [0,1] defined 

recursively by 

50 ( 0 ;) = and Sk+i{x) = Tk^i(Sk(x)) 

where acts only of the fc-th digit: Tk{x) = x — (TV — 1) ■ N~^ if rc G Cj- and Tk{x) = x N~^ 
if a: G [0,1] \ Ck- The sets Ck are defined by 


c. - y {jiL 


a=i 


2 _ L_ 

J^k ’ J\[k-1 


Moreover, it is easy to see that 


(T/c){{>C[o,i] = '^[0,1] ^ Sk S uniformly. 


Hence 

j f{x)SiC'^ = J f(x)dx, V/ e Co([o, 1]). 

Now, it remains to show that T is optimal. But this is true thanks to the fact that \A.2\ implies 
that the plan induced by T satisfies the hypothesis of Lemma |4.3| 

It is clear that we can reduce to prove the non-differentiability and the fractal properties only for 
the map S. The non-differentiability comes from the fact that for each 2 : G [0,1], we have that the 
base- TV representation 2 = is such that can’t be definitively TV— 1. Now it is sufficient 

to choose those kj such that ak^ 7 ^ TV — 1 and consider the numbers Zj = z (TV — 1 — afcj )/TV^J 
and z'^ = z^ 1/TV^-? (depending on whether Okj = 0 or not) which have the same digits as 2 apart 
from the digit. Then it is straightforward to see that S{zj) — S{z) = —(a^^ -h 1)/TV^J while 

S{z'2 — S{z) = z'^ — 2 ; in particular we have ^ ^ = 1 while ^*^5 

letting j ^ 00 , we get that S is not differentiable at the point 2. 

As for the fractal property: we consider the transformation 


^(9)(x) 


jfg{Nx - i) + ^ 

^g(Ar^_7V_l) 


for ^ < X < , with i = 0,..., N — 2 

for < X < 1. 


In order to see what the construction is doing we imagine to divide [0,1]^ in a grid X of 
squares and then we are putting scaled copies (by a factor TV) of the original function in the above 
diagonal squares and in the rightmost bottom one. Now it is easy to see that S is the uniqu(0 
fixed point of J^, and this gives the property of self-similarity. □ 


Remark 4.7. Another proof of Theorem |4.6| can be done noticing that the transformation J-^ 
(the one defined in the proof) acts as a 1/TV-contraction in the space of measure preserving L°°- 
bijections from [0,1) to [0,1). 

Corollary 4.8. Let /r = uniform measure on the d-dimensional cube [0,1]'^ in 

R^ and suppose c : ([0,1]^)-^ R the N-dimensional repulsive harmonic cost 

N N 

c(xi,... ,xjv) = - \xj-Xi\‘^, (a;i,...,Xiv) e ([0,1]'*)^ 

i = l j = 

Then, there exits an optimal cyclical transport map T : [0,1]*^ ^ [0,1]^; in particular 

min f cd'y = min [ c(x,T{x),T^^^ (x),... ,T^^~^\x))dp. 

7Gn^(/i)7 y 

J-C-N) _J 

Moreover, T is not dijferentiable at any point. 


Proof. As we already observed, the problem with the cost c is equi valent to the problem with the 
cost \xi -h ... + xjsfl^ and the result follows from the Theorem |4.6[ □ 


^We notice that .J" is a contraction: ||.!F(g') — .F’(5^)j|oo ^ ~ sdloo 
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We remark the construction of T in the proof of the Theore m |4.6| also works if TV = 2 and T 
is exactly the optimal transport map described in the example l |4.2[ |. 

The unexpected aspect of Corollary |4.8| is the existence - for TV > 3 - of an optimal transport 
map which is not differentiable almost everywhere. It turns out that this optimal map T could 
be not unique if d > 1 and TV > 3 and, in that case, we can construct explicitly a regular optimal 
map. This kind of richness appears already in the case d = 1 with the Coulomb cost; however in 
that case we have uniqueness if we restrict ourselves to the symmetric optimal plan. Here this is 
not the case as it is easy to see that modifying the action S on the digits (the important thing is 
that when we see tS as a permutation, it is a cycle), we obtain another map, and the symmetrized 
plan is not equal to the one generated by the map described before. 

In the following, we are going to present some concrete examples where we have other explicit 
solutions for such kind of optimal maps for (AifCweak) in ( |4.1[ |. For these examples the goal is to 
show that there can be smooth optimal maps but they don’t necessarily satisfy the “group rule”, 
that is, we can find maps T 2 ,...Tjv such that (/d, Ti,T 2 ,..., Tjv—i)jipi is an optimal plan but 
Ti o Ti 7 ^ T-i for any index i G {1,..., TV — 1}. 

Example 4.9 (3 particles, asymmetric). Consider the case when three particles are distributed in 
as Gaussians, fii = p ,2 = ^ 3/2 +^2 + ^ 3 )) “ 2 +^2 + ^ 3 ))- 

In this case, we can verify that the couple (T, 5) of maps T,S : R*^ R^, T{x) = x a.e. and 

S{x) = —2x a.e. is admissible and it is optimal since x + T{x) + S{x) = 0. 

Example 4.10 (2N particles on 5^, ). Suppose /ii,..., p 2 N uniform probability measures on the 
circle . The rotation map Rq : R^ ^ R^ with angle 0 = tt/TV is an optimal transport map. 
Also, the maps Rke, fc = 2,..., TV, are optimal transport maps for the repulsive harmonic cost. 

Example 4.11 (2N particles, breathing map). Suppose pi,..., p 2 N uniform probability measures 
on C R^. Consider the vector u, A : R^ —>■ R^ the antipodal map A{x) = —x and Rg : R^ —>■ R^ 
the rotation of angle 9 = tt/TV and direction v. Then, T = Ao Rg is a cyclic optimal transport 
map. 

The optimal solution 7 = (x, T, ,..., )jjpi is called “breathing” solution |83| : the 

coupling 7 represent the configuration where the 2TV electrons are always at the same distance 
from the center, opposite to each other in the equilibrium configuration. Notice that, the map 
G = is also an cyclical optimal transportation map. Moreover, T is a function and it is a 
gradient of the convex-concave function (l){xi,X 2 ,..., ^ 2 iV— 1 ? ^ 2 N^) = 5(^1 ‘ ‘ ‘ ~ 

^{x^ + xl H-l-a;2^). 

Example 4.12 (2N particles in R^, symmetric p). In this case, on can consider the maps T{x) = x 
a.e. and S{x) = —x a.e. which are such that 

c{T{x), S{x),... ,T{x), S{x)) = 0. 

We notice that in this case we have (t) = x and so in particular this solution has the cyclic 
structure 

iT{x), S{x),T{x), S{x)) = (x, Six), ,.. _ 

Example 4.13 (Optimal Maps which doesn’t satisfy a group law). Let us consider the case d = 1 
and N = 3 where the measures are pi = p ,2 = = ^^|[—i,i]- 

Then we define the maps 

. (X1 if a: < 0 . f —1 — 23 ; if X < 0 

T{x) = { - S(x) = { 

1^3: —1 if 3 : > 0, 1^1 —2t if X > 0. 

We have that Typ = p and S^p = p, and moreover x -\-T{x) S(x) = 0; in particular {Id, S, T)^p 
is a fiat optimal plan and so the thesis. 


Exam ple 4.14 (Optimal Transport diffuse plan). Let us consider the same problem as in Example 
4.13 Now we consider a general symmetric plan 7 = f{f3iax{\x\,\y\,\z\}), where H is 

defined as H = {x y + z = 0} D {| 3 :| < 1 , \y\ < 1 , \z\ < 1 } and / is a function to be chosen later. 
This is a symmetric flat optimal plan; now we compute the marginals. Since it is clear that 7 is 
invariant under x 1 —>• —x, it is sufficient to consider the marginal on the set 3 : > 0. But then we 
can make the computation 
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Jo Jo 


V 30 (x)/(max{|x|, |;/|, \x + y\}) dxdy 

i>0 ~ 

1 — x rl rO 

y/^(l){x) f {x-\-y) dy dx-\-y/Z / / (^{x) f {x) dy dx 

Jo J-X 

(l>{x)f{\y\)dydx 


■■ j V 3 (f){x) ^xf{x) + 2 j f{t) dt^ 


dx. 


In particular the choice f{x) = ^x gives the marginals equal to We notice also that 

any even density p = h{\x\) for some decreasing function h : [0,1] [0, oo) can be represented in 

this way: in fact it is sufficient to choose f{x) = ~ fx • 

Example 4.15 (A Counterexample on the uniqueness > 3). As mentioned in theorem |2.16| 
and in m, it was already understood by Pass that in these high dimensional cases, the solution 
of repulsive costs may also be non unique, as opposite to the two marginals case. On a higher 
dimensional surface there can be enough wiggle room to construct more than one measure with 
common marginals, as shown in the examples |4.9|{4!T^ In most of the cases, the non-uniqueness 
seems to be given by the symmetries of the problem, but in Example |4.14| and Corollary |4.8| this is 
not the case, as we exploit the fact that the dimension of the set c{x^y-\-z) — (l){x) — (p{y) — (l>{z) = 0, 
is greater than the minimal one. 

Finally, the last proposition of this section states that when d = 1 and odd N we have no hope 
in general to find piecewise regular cyclic optimal transport maps. 

Proposition 4.16. Let y = >C|[o,i] N > 3 be an odd number. Then, the infimum 

c{x,T{x),T(^\x),...,n^-^^)dti : j. 

is not attained by a map T which is differentiable almost everywhere. 

Proof. First of all we notice that if T is differentiable almost everywhere then also has the 
same property, thanks to the fact that T^y = y. In particular, since the Lusin property holds true 
for for every i = 1,..., N and T^y = y, the change of variable formula holds and in particular 
we have that (T^*^)'(x) = ±1 for almost every x (notice also that T is bijective almost everywhere 
since T^^'^{x) = x). 

Since y is A’-flat we have that the condition on T in order to be an optimal cyclical map is 
X -h T{x) + (x) + ... + (x) = A/2; now we can differentiate this identity and so we will 

get 

1 + T'{x) + (T(2))'(x) + ... + = 0 for a.e. x. 

But this is absurd since on the left hand side we have an odd number of ±1 and their sum will 
be an odd number. □ 

In conclusion, also in the case of the repulsive harmonic cost, the picture is far from being 
clear: an interesting structure appears when {yi}^i is a fiat A-tuple of measures but we still 
can’t characterize this property. Moreover, in the fiat case in which yi = y, for example when 
y = we have both diffuse optimal plan and a cyclical optimal map. 

An interesting open problem is whether for any A-flat measure y, say absolutely continuous 
with respect to the Lebesgue measure, we have a cyclical optimal map and a diffuse plan. 


5. Multi-marginal OT for the Determinant 

We are going to give a short overview of the main results in El, where Carlier and Nazaret 
consider the following optimal transport problems for the determinant: 


(MlCnet) sup 

76 n{{K‘^)'^ 


p / 


det(xi,..., Xa)d'y{xi, ...,xa) 


( 5 . 1 ) 
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and 

{MK.\oe.t\) sup / |det(Ti,...,3:d)M7(^i>-*-,a:d), (5.2) 

where fii,..., fid are the uniform probability measures in In addition, in order to guarantee 
existence of a solution, we assume that there exist pi,... , E [ 1 , oo[ such that 


^7 = ^’ 


and 


d p 

5 /. 


\xi\^ 


H'i Pi 


-d^{xi) < + 0 O. 


Notice that for this particular cost the problem makes sense only when N = d. In 

the following, we will focus on proble m (|5.1[ ) and exhibit explicit minimizers 7 in the radial 
case. Clearly, the difference between | |5.1| l and ( |5.2[ l is that the second one admits positively 
and negatively oriented basis of vectors, while the first one “chooses” only the positive ones. 
Moreover, if we assume that among marginals fii,..., fid there exist two symmetric probability 
measures fii^fij^i 7 ^ j, i.e. fii = (-Id)tj/^i and fij = (-Id)[j/j,j, then any solution 7 of satisfies 

det( 3 :i,... ^Xd) > 0 7 -almost everywhere and so solves also ( |5.2[ l (Proposition 6 , 1171 ). 

Similarly to the Gangbo-Swi§ch cost the Monge-Kantorovich problem for the determinant 
can be seen as a natural extension of classical optimal transport problem with two marginals 
and so, it is equivalent to the 2-marginals repulsive harmonic cost ( |4.1[ |. 

Indeed, we can write in the two marginals case, det(xi,a: 2 ) = (xi,Rx 2 ), where 7? : R^ ^ R^ is 
the rotation of angle —7rl2. Hence, since fii and fi 2 have finite second moments, up to a change of 
variable X 2 = Rx 2 , the problem {AiKoet) iu is equivalent to the classical Brenier’s optimal 

transportation problem: 


argmax / det{xi,X 2 )d'y{xi^X 2 ) 


= argmax 

7en(Aii,/i2) 

= argmax 

7en(Aii,A2) 

= argmin 

7en(^i ,ii2) 


/ {xi,X2)d'y{xi,X2) 

/ {x\,X2)d^ 

Jb? 

L 


In - X2\^ 


d'y{xi,X2) - C 


where C = l/2(/ l^ip dpi + f \x 2 \‘^ dp 2 ) and p 2 = RftP 2 - 

In the sequel, we are going to construct maximizers for (|5.1| t, thanks to some properties of the 
Kantorovich potentials of the dual problem associated to ( |5.1| l (see theorem |5 .1 1 bellow), 


(ICT) inf 


Ui{xi)dfii{xi) : det(Ti, 


,^d) < ^Ui{Xi) 


(5.3) 


In |17| . the authors provide a useful characterization of optimal transport plans through the 
potentials Ui, given by theorem |5.1| . In addition, by means of a standard convexification trick we 
obtain regularity results on the Kantorovich potentials. 


Theorem 5.1. A coupling 7 E n((R'^)'^, /j,i,..., fid) is optimal in if and only if there exists 

lower semi-continuous convex functions Ui : R^ ^ R U { 00 } such that for all i E {1,..., d}, 


d d 

Axj), on (r7; 
i=l i=l i^j 

d d 

Ui{xi) = ^^u*((— 1 )^"*"^ Xj), 7 — almost everywhere; 

i=l i=l 

(— 1 )^"^^ f\xj ^ dui{xi), 7 — almost everywhere. 

where A?=i ‘^3 denotes the wedge product and, for every i, u* is the convex dual of the Kantorovich 
potential Ui . 
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Now, the main idea is to use the geometrical constraints on the Kantorovich potentials Ui 1^, 
given by the theorem in order to construct an explicit solution. 

We illustrate the theorem (53 and explain how to construct a particular optimal 7 for (H^by 
an example in the three marginals case. Let = piC^, i = 1 , 2 , 3 radially symmetric probability 
measures on the 3 -dimensional ball B. 

In this particular situation, the optimizers of and have a natural geometric inter¬ 

pretation: what is the best way to place three random vectors x,y^z, distributed by probability 
measures pi, P 2 -, Pz 011 fhe sphere, such that the simplex generated by those three vectors (x, y, z) 
has maximum average volume? 

Suppose 7 G I\.{B, pi, p2, pz) optimal in when d = S. Prom optimality of 7, we have 

ui (x) -h U2{y) + us{z) = det{x, y,z), 7 — almost everywhere. 

Applying the theorem we get 

( 'u,2{y) + U2,{z) = u\{y/\ z) 

< ui{x) + uz{z) = U2{—x A 2) , 7 — almost everywhere, 

I Ui (x) U2 {y) = u^{x A y) 

and, 

{ Vui(x) = y A z 

Vu2(t/) = —X A z , 7 — almost everywhere. ( 5 - 4 ) 

\Juz{z) = X Ay 

It follows from ( | 5 . 4 [ |, given a vector x on in the ball, the conditional probability of y given x is 
supported in a “meridian” M{x) 

M{x) = {y e 5 ^ : {SJux{x),y) = 0 }, 

where 5 ^ is the 2 —sphere. Finally, assuming that (x, Vui(rr)} 7^ 0 , the conditional probability of 
2 given the pair {x, y) is simply given by a delta function on 2 

Vui{x) A Vu 2 {y) 

{x, Vui(x)} 


In particular, we have 

(x, Vui(a;)) -h {y,Vu2{y)) -h (2, Vu3(2)} = det{x, y,z) = det('Vui{x),'Vu2{y),Vu3{z)). 

Example 5.2 (An explicit solution in the ball B C R^). Suppose pi = C^, i = 1 , 2 , 3 , the 
3 -dimensional Lebesgue measure in the ball B GR.^. The following coupling 7* 

/ = ralm f ( f /(a:, ^^^ ']dx, V/g C(B^,R). ( 5 . 5 ) 

Jb 3 C 3 {B) Jb \Jm(x) 27 r J 

is an optimizer for with d = 3 . Indeed, from we can show explicit potentials u^(x) = u^ix) = 
u^{x) = |a:p/ 3 ; clearly we have 

det{x,y,z) < |a:||y||2| < |3;|®/3 + |y |®/3 + \z\^/‘i, V {x,y,z) e B, 

with equality when |x| = \y\ = \z\ and x, y, z are orthogonal. Since 7* is concentrated on this kind 
of triples of vector we have the optimality. Finally, by a suitable change of variable, it is easy to 
see that 7* G IVz{E\). 

Some comments on the radially symmetric d-marginals case: In ca. for d radially symmetric 
probability measures, the authors exhibit explicit optimal couplings 7*. In their proof, two aspects 
were crucial: the first one is remark that if {pi}f—i are radially symmetric measures in then 
the optimal Kantorovich potentials Ui(xi) = are also radially symmetric; in particular the 

system ( | 5 . 4 [ | for general d implies that the support of 7* is contained in the set of an orthogonal 
basis. The second observation is to notice that in the support of 7* we have Hi{\xi\) = \xi\, where 
Hi is the unique monotone increasing map such that (H*)^pi = pi, where H*(x) = |^iLi(|a:|). 
This is done analyzing the correspondent radial problem (with cost c(ri,..., r^) = n • • • r^), using 
the optimality condition (pi{ri) = diC and the fact that in this case ri(j)'^{ri) = ri0'^(ri) = c > 0; 
then, exploiting the convexity of we get ri = Hi{ri) for some increasing function Hi, that is 
uniquely determined. 





OPTIMAL TRANSPORTATION THEORY WITH REPULSIVE COSTS 


27 


Existence of Monge-type solutions: In the 3 -marginals case in the unit ball, by construction of 
the the coupling 7* in or, more generally, the optimal coupling in the d-marginal case (see 

theorem 4 in [ 17 ) 1 . we can see that their support are not concentrated in the graph of cyclic maps 
T, or simply on the graph of maps Ti,... as we could expect from corollary 

| | 2 . 15 [ ). In other words, 7* in is not Monge-type solution. 

The existence of Monge type solutions for the determinant cost is still an open problem for 
odd number of marginals. From the geometric conditions we discussed above, in the case in which 
pi = fi a radial measure, if Monge solutions exists then, for every x E (a^, Ti (x),..., Tci_i{x)) 
should be an orthogonal basis, and |Ti(x)| = \x\ , z = 1 ,. .., d — 1. 

For the interesting even dimensional case, we can observe a similar phenomena remarked in 
the repulsive harmonic costs, concerning the existence of trivial even dimensional solutions for the 
Monge problem in We expect the existence of non-regular optimal transport map 

also to this case. 

Example 5.3 (Carlier & Nazaret, [E])- The even dimensional phenomena: As in the repulsive 
harmonic cost, it is easy to construct Monge minimizers for the determinant cost for even number 
of marginals > 4 . For instance, suppose c(xi, X2,2:3,0:4) = det{xi^X2,xz,X4)y define transport 
maps Ti, T2, T3 : 5 ^ by, for x = (a;i, 3^21 ^3? ^^4) £ B 


/ -X2 \ 




/ -X4 \ 

Xi 

, T2{x) = 

X4 

. T3(x) = 

-X3 

— 3:4 


3:1 


X2 

\ ^3 / 


\ -^2 / 


\ XI J 


We can see that 77^ = ( 7 d,Ti,T2,T3)(t7 is a Monge-type optimal transport plan for and 


6. Numerics 

Numerics for the multi-marginal problems have so far not been extensively developed. Dis¬ 
cretizing the multi-marginal problem leads to a linear program where the number of constraints 
grows exponentially in the number of marginals. In m Carlier, Oberman and Oudet studied the 
matching of teams problem and they were able to reformulate the problem as a linear program 
whose number of constraints grows only linearly in the number of marginals. More rencently a 
numerical method based on an entropic regularization has been developed and it has been applied 
to various optimal transport problem in 0 1101129| . Let us mention that the Coulomb cost has 
been treated numerically in various paper as in [Mills], where the solutions are based on the 
analytical form of transport maps given in ESI , in EH, where parameterized functional forms of 
the Kantorovich potential are used, in |201119| . where a linear programming approach has been 
developed to solve the 1 -dimensional problem, and, as already mentioned, in Ho]. In this section 
we focus on the regularized method proposed in [^ 110112 ^ and we, finally, present some numerical 
experiments for the costs studied above. 

6.1. The regularized problem and the Iterative Proportional Fitting Procedure. Let 

us consider the problem 

{A4K,) inf / c{xi,... ... ,xjs[)dxi,... ,dxN, ( 6 - 1 ) 

where N is the number of marginals pi which are probability distributions over R*^, c{xi, ...,xn) 
is the cost function, 7 the coupling, is the probability distribution over R^^ 

Remark 6.1. From now on the marginals pi and the coupling 7 are densities and when the optimal 
coupling 7 is induced by maps T^, we write it as 7 = pi{xi)5{x2 — T 2 {x\)) • • • 5{x]s[ — Tn{x\)). 

In order to discretize ( |6.1| |, we use a discretisation with points of the support of the fcth 
marginal as ■ ,M^- If the densities p}^ are approximated by p^ ~ 

min jjv- ( 6 - 2 ) 

where 11^ is the discretization of 11, = c{xj^,--- ,Xjj^) and the coupling support for 

each coordinate is restricted to the points ,Mrf thus becoming a matrix again 
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denoted 7 with elements The marginal constraints Ci (such that H/c = f]^—iCk) 

becomes 


Ck ■■= { 


7 e 


= I,-- - , A^d}- 


E 

As in the continous framework the problem ( | 6 . 1 [ | admits a dual formulation 

N M 


( 6 . 3 ) 


max 

A; = lj = l 
N 

S.t. Ujj^ < 

k = l 


( 6 . 4 ) 


V ifc = 1 , • • • , Md, 


where Uj^ = Uk(xjf^) is the kth Kantorovich potential. One can notice that the primal l | 6 . 1 [ l 
has {Md)^ unknown and Md x N linear constraints and the dual problem ( 6 . 41 has Md x N 
unknown, but (Md)^ constraints. This actually makes the problems computationally unsolvable 
with standard linear programming methods even for small cases. 


Remark 6 . 2 . We underline that in many applications we have presented (as in DFT) the marginals 
are equal (/i.j^ = fij, V k £ { 1 , • • • ,N}). Thus the dual problem can be re-written in a more 
convenient way (but still computationally unfeasible for many marginals) 


M 

max E Nujtij 

" 3=1 

N 

S.t. 55 '^Jk — '^3k = ^ 1 ’ ’ ’ ) ^di 

k = \ 


( 6 . 5 ) 


where Uj = Uj^ = u{xj^). Now the dual problem has Md unknown, but {Md)^ linear constraints. 

A different approach consists in computing the problem ( | 6 . 1 [ l regularized by the entropy of the 
joint coupling. This regularization dates to E. Schrodinger I8II and, as mentioned above, it has 
been recently introduced in many applications involving optimal transport 1101 [ 2 ^ I 40 | . Thus, 
we consider the following discrete regularized problem 

E +^^(0-) (6-6) 

7fcC . 


where E('y) is defined as follows 

£(7) = if7>0 

1+00 otherwise, 

and C is the intersection of the set associated to the marginal constraints (we remark that the 
entropy is a penalization of the non-negative contraint on 7). After elementary computations, we 
can re-write the problem as 

minA'L(7|7) (6.8) 

_^ 7 i i AT 

where KL{'y\^) = ^ 7ii,...7^ log(3——) is the Kullback-Leibler distance and 

^ 

_ Jn 

= e e . ( 6 . 9 ) 

As explained in section when the transport plan 7 is concentrated on the graph of a trans¬ 
port map which solves the Monge problem, after discretisation of the densities, this property is 
lost along but we still expect the matrix 7 to be sparse. The entropic regularization spreads the 
support and this helps to stabilize the computation as it defines a strongly convex program with 
a unique solution 7^. Moreover the solution 7^ can be obtained through elementary operations. 
The regularized solutions 7^ then converge to 7* (see figure]^, the solution of l | 6 . 1 [ l with mini¬ 
mal entropy, as e ^ 0 (see m for a detailed asymptotic analysis and the proof of exponential 
convergence). 
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Figure 2 . Support of the coupling 7^ for the Coulomb cost and pi = = 

(1 + cos(7rx/2))/2. The simulation has been performed on a discretization of 
[— 2 , 2 ] with Md = 1000 . 

In order to introduce the Iterative Proportional Fitting Procedure (IPFP), we consider the two 
marginals problem 

minY^Cij7ij. (6.10) 

7GC 

The aim of the IPFP is to find the KL projection of 7 on the set C = {'fij S X R*^<i : 

u, 7j^ = n {'yij G X “ ^j}- writing down the lagrangian associated 

to | | 6 . 10 | > and computing the optimality condition, we find that 7 *^ can be written as 

7 j^ = aibj^ij with ai = bj = ( 6 - 11 ) 


where ui and Vj are the regularized Kantorovich potential. Then, Uj and bj can be uniquely 
determined by the marginal constraint 


ai = 




E, bj'Tij '' 


fo. - _7_ 

^3 ~ - • 

<^ai3 


Thus, we can now define the following iterative method 


Ln+l _ ^3 n + l _ 


( 6 . 12 ) 


(6.13) 


Remark 6.3. In Ea Ruschendorf proves that the iterative method ( |6.13| | converges to the KL- 
projection of 7 on C. 


Remark 6.4. Ruschendorf and Thomsen (see [78]) proved, in the continous measure framework, 
that a unique KL-projection exists and takes the form 'y{x,y) = a(x) 0 b{y)^{x,y) (where a{x) 
and b(y) are non-negative functions). 


The extension to the multi-marginal framework is straightforward but cumbersone to write. It 
leads to a problem set on N M^f-dimensional vectors ajp^ ^, j = 1, ■■■, = 1, ■ ■ ■ , ^d- 

Each update takes the form 




Pij 


y. 






a!*+^ a:i+K..a"+} . a" , . ...a" . 

1,^1 2,^2 j — + 1 
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Example 6.5 (IPFP and 3 marginals). In order to clarify the extension of the IPFP to the multi¬ 
marginal case , we consider 3 marginals and we write down the updates of the algorithm 


n+l 

Pil 


l,ll 

^2,22 

Uo ■ 
3,t3 

- 

Pi2 


2,^2 

Eii,i3 7'ii72;*3 

Uo - 

3,^3 

n + l 

Pi3 


3,23 

- n-|-l 

2^iiP2 7nd2d3 “i,2i 

aO+i 

2,^2 


In the following sections we present some numerical results obtained by using the IPFP. 

6.2. Numerical experiments: Coulomb cost. We present now some results for the multi¬ 
marginal problem with Coulomb cost in the real line. We consider the case where the N marginals 
are equal to a density p (we work in a DFT framework so the marginals rapresent the electrons 
which are indistinguishable). We recall that if we split p into N pi with equal mass {f pi(x)dx = 
f p(x)dx), then we expect an optimal plan induced by a cyclical map such that T^pi = Pi+i 
2 = 1 , • • • , — 1 and T^pj;^ = pi. 

The simulations in figure[^are all performed on a discretization of [—5, 5] with = 200, with 
marginals pi = p(x) — ^(1 + c.os{^x)) i = 1, • ■ ■ , and e = 0.02. If we focus on the support 
of the coupling 712 ( 2^5 y) = (^i 5 ^ 2 ) 7 ( 3 ^! y?- 2 ) we can notice that the numerical solution correctly 
reproduces the prescribed behavior: the transport plan is induced by a cyclical optimal map (see 
section]^. We refer the reader to [TO] and m for examples in higher dimension. 

Remark 6.6. Theorem |3.7| actually works also for other costs functions as c(xi,--- ,xn) = 
—log{\xi — Xj\). This means that if we use the density p{x) = ^(1 cos{^x)), we expect 

to obtain the same solution as for the Coulomb cost. Thus if we now consider the 3-marginals 
case and the discretized cost Cji,j 2 ,j 3 — ~ ^^ 2 !) ~ ~ ^jsl) ~ ^^ 9 i\^j 2 ~ 

we can notice (see figure]^ that we recover the same optimal coupling as the corresponding case 
for Coulomb (first row of figure]^. The simulation is performed on a discretizetion [—5,5] with 
Md = 200, with marginals p{x) = ^(1 + cos{^x)) N = 3 and e = 0.01. 


6.3. Numerical experiments: repulsive Harmonic cost. As shown in section the min- 
imizers of the OT problem with the harmonic cost are als o m inimizers of the problem with 
c(a:i, • • • ,Xjv) = 1^1 + • • • + so for the discrete problem ( 6.2 1 we take = \xj^ + • • • + 

r p. Let us first consider the two marginal case and the uniform density on [0,1] (see example 
and, as expected, we find a deterministic coupling given by 7 (x, y) = p{x)5(y — T(x)) where 
T{x) = 1 — x, see figure]^ The simulation in figure]^ has been performed on a discretization of 
[0,1] with Md = 1000 gridpoints and e = 0.005. 

The multi-marginals case is more delicate to treat: the original OT problem does not admit 
a unique solution whereas the regularized problem does. As we have explained in section |6.1| 
the regularized problem is a strongly convex problem which admits a unique solution and the 
resulting coupling is the one with the minimal entropy. However we are able to make the IPFP 
algorithm converge to a select ed coupling among the optimal ones for the original problem. Let 
us focus on the example 4.13 In this case we have all marginals equal to p = and we 

can find a deterministic coupling given by 7 ( 2 :, y, 2 ) = p{x)5{y — T{x))5{z — S{x)) for the maps 
T{x) and S{x)). In order to select this coupling, the idea is to modify the cost function by adding 
a penalization term p( 37 , y) = which makes pi and p 2 be as far as possible (if we consider 

pi and p 2 as particles). Thus, the discretized cost now reads as 


= \xj^ -h: 




(6.14) 


kil - ^j2 I 

Then we expect, as shown in figure]^ that the projection 712 ( 3 ;, j/) = (ei, 62 ) 7 ( 3 ;, j/, 2 :) and 
713 ( 0 ^, 2 :) = (ei, 63 ) 7 ( 0 :, y, 2 ) of the computed coupling are induced by map T and S respectively. 
The simulation has been performed on a discretization of [—1,1] with Md = 1000 gridpoints, 
e = 0.0005 and y = 0.1. 

We, finally, take the same marginals as in the previous example, but we do not add the penal¬ 
ization. In this case we expect a diffuse plan since, if we forget about the marginal constraint, it 
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Figures. OntheLeft: support of the optimal coupling'j 12 {x,y), on the Right: 
graph of the optimal coupling •yi2{x,y) , in the cases = 3 ( 1 st Row), = 4 
( 2 nd Row), N = b ( 3 rd Row). The dashed lines delimit the supports of jli , 
with 2 = 1, • • • , N. 



Figure 4 . (logarithmic cost) Left: support of the optimal coupling •yi2{x,y) 
for N = S. Right: support of the optimal coupling 712(^5 for N = 3 . The 
dashed lines delimit the intervals where jli , with 2 = 1, • • • , 3 , are defined. 
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Figure 5. Left: support of the optimal coupling. Right: optimal coupling. 



Figure 6. Top-Left: support of the optimal coupling ^i2(x,y). Top-Right: 
optimal coupling ^i2{x,y). Bottom-Left: support of the optimal coupling 
7i3(x,2:). Bottom-Right: optimal coupling ^is{x, z). 
is simple to show that the second order F-limit of the regularized problem as e ^ 0 is the relative 
entropy of 7 with respect to 'H'^\{x-\-y-\-z=o}i expect that the F-limit is the same also when we 
add the marginal constraint. 

In figure]^ we present the projection 712(^5^) = (^i, 62)7(3:, y, 2:) of the computed coupling 
(it is enough to visualize only this projection because of the symmetries of the problem) and, as 
expected, we observe a diffuse plan. . The simulation has been performed on a discretization of 
[— 1 , 1 ] with Md = 1000 gridpoints, e = 0 . 0005 . 


6.4. Numerical experiments: Determinant cost. Numerical simulations for the multi-marginal 
problem with the determinant cost present an obvious computational difficulty: in order to com¬ 
pute the solution just for the 3 -marginals case, we have to introduce a discretization of (as 
the cost function now is c(a:i,--- ,xn) = det{xi, • • • ,X]\j) with Xi £ R^). However if we take 
radially symmetric densities pi as marginals (see section]^ for a more complete description of the 
problem in the measure framework), Carlier and Nazaret show that (AiK.) can be reduced to a 
1-dimensional problem: the only unknowns are the relations between the norms ri = ||3:i|| of each 
vector. These relations can be obtained by solving the following problem 

iV 

min - / (FT n)7(ri,--- ,riv)dri ■••dr^, ( 6 . 15 ) 

,Xn) J f-T 
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Figure 7 . Left: support of the optimal coupling 712(3;, y). Right: optimal 
coupling 712 ( 3 ^! y)- 

where Ai(r) = a{r)pi(T) (e.g. TV = 3 then r = {r^9^(p) and cr(r) = 47rr^). Moreover, for this 
problem we know that there exists a unique (deterministic) optimal coupling (see Proposition 5 
in El). The discretized cost now reads consider the 3-marginals 

case and all densities equal to the uniform density on the ball So.5 = { 3 ; G M^|||3;|| < 0.5}. As 
proved in [13, the optimal coupling is actually supported by the graph of maps which rearrange 
measure (see figure]^). The simulation in figure has been performed on a discretization of 
[0, 0.5] with Md = 100 and e = 0.01. In the same way we can take all marginals Aj(r) = 47rr^e“^’' 



Figure 8 . (Uniform density) Left: support of the optimal coupling 712 ('^i, '^2)- 
Right: optimal coupling ^ 12(^1 ^r 2 )- 


and we obtain again a deterministic coupling, fij 
discretization of [0, 3] with = 300 and e = O.l 



Figure 9 . (Exponential density)Left 
7 i 2 (?"i 5 ?" 2 )- Right: optimal coupling 


The simulation has been performed on a 



support of the optimal coupling 
iri,r2). 


We, finally, present a simulation (see figure |ll[ l with three different marginals Xi (see figure 
|10[ l. As one can observe we obtain that the projection of the optimal coupling are concentrated 
on the graph of a map which rearranges the densities in a monotone way. The simulation has 
been performed on a discretization of [0,4] with = 300 and e = 0.07. 
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Figure 10. Densities Xi, X 2 and A3. 



Figure 11. Top-Left: optimal coupling ■y 12 {ri,r 2 ). Top-Center: optimal cou¬ 
pling ^i 2 ,(ri^r^). Top-Right: optimal coupling ^ 2 ziT 2 -,'f'z)- Bottom: supports 
of those couplings. 

7. Conclusion 

In this survey we presented a general famework for the multi-marginal optimal transportation 
with repulsive cost. We also gave some new results which highlight the importance of the multi¬ 
marginal approach. As mentioned there are still a lot of open questions about the solutions for 
the problems presented and we hope to treat them in future works. 
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